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The superfluid/normal-fluid interface of liquid "^He is investigated in gravity on earth where a small 
heat current Q flows vertically upward or downward. We present a local space- and time-dependent 
renormalization-group (RG) calculation based on model F which describes the dynamic critical ef- 
fects for temperatures T near the superfluid transition Tx- The model- _F equations are rewritten 
in a dimensionless renormalized form and solved numerically as partial differential equations. Per- 
turbative corrections are included for the spatially inhomogeneous system within a self-consistent 
one-loop approximation. The RG flow parameter is determined locally as a function of space and 
time by a constraint equation which is solved by a Newton iteration. As a result we obtain the 
temperature proflle of the interface. Furthermore we calculate the average order parameter {ijj), 
the correlation length ^, the specific heat Cq and the thermal resistivity px where we observe a 
rounding of the critical singularity by the gravity and the heat current. We compare the thermal 
resistivity with an experiment and find good qualitative agreement. Moreover we discuss our previ- 
ous approach for larger heat currents and the self-organized critical state and show that our theory 
agrees with recent experiments in this latter regime. 
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I. INTRODUCTION 

On earth in liquid ^He the gravity is an external force 
which causes a space dependent pressure p = p{z) de- 
pending on the altitude coordinate z. Since the critical 
temperature of the superfluid transition T\ = T\{p) de- 
pends on the pressure p, in the helium the critical tem- 
perature T\{z) ~ T\{p{z)) varies with the altitude z. In 
leading approximation it is a linear function of the alti- 
tude 



niz) ^ Tx{zo) + (dTx/dz) {z - zo) 



(1.1) 



where the gradient is determined experimentally as [ij 
dTx/dz ~ +1.273 /xK/cm. The sign is positive which 
means that the critical temperature increases with the 
altitude z. 

In thermal equilibrium, the local temperature of the 
helium r(r, t) =T is constant with respect to any space 
and time variable. If in an experiment we choose the 
temperature T = T\{zq) wc find an interface at z = zq, 
which separates superfluid ^Hc in the upper region z > 
zo where T < T\{z) from normal-fluid '*IIc in the lower 
region z < zq where T > T\{z). This interface is the 
main concern of the present paper. 

Correlation effects imply an interface which is not 
sharply defined but smeared out over a certain length 
scale ^g. Ginzburg and Sobyanin [2[ have calculated the 
order parameter profile V'(^) for hquid ^He in gravity 
within their theory which is a mean-field theory mod- 
ified by scaling functions in order to incorporate the ef- 
fects of critical fiuctuations and the critical exponents to 
some extent. They find the characteristic length scale 
ig = 67 /^m (see Fig. 4 and Eq. (3.49) in Ref. i). 

A heat current Q flowing from bottom to top in 
the direction of z enhances the formation of the 



superfluid/normal-fluid interface. Heat transport phe- 
nomena imply a space-dependent temperature T(z) with 
a negative gradient dT/dz < which acts opposed to 
the positive gradient of the critical temperature Ta(z). 
Onuki d, 0] has investigated the interface under a heat 
flow Q within a dynamic mean-field theory modified by 
scaling functions. He finds that the thickness of the inter- 
face decreases according to ^ Q~^l^ with increasing 
heat current Q. 

While on earth the gravity acceleration g = 9.81 m/s^ 
is constant; the heat current Q can be varied in the ex- 
periment. For large heat currents Q > Qq, the heat- 
current effects dominate, where on the other hand for 
small heat currents Q < Qo gravity effects dominate. 
The heat current which separates both regimes, is about 
Qo = 70nW/cm^. In this paper we focus on small heat 
currents where gravity is the main effect. 

The critical dynamics of liquid *He near the superfluid 
transition T\ is described by a hydrodynamic model with 
Gaussian fluctuating forces which is called model F in the 
classification of Hohenberg and Halperin [5| . This model 
has originally been derived by Halperin, Hohenberg, and 
Siggia [6j in order to describe the critical dynamics of 
a planar ferromagnet, which is in the same universality 
class as liquid ^He. The field-theoretic renormalization- 
group theory of model F has been elaborated by Dohm 
0, The specific heat and the thermal conductivity 
have been calculated up to two-loop order and com- 
pared with very accurate experimental data 0, [l3| • In 
this way, the renormalized coupling parameters and some 
other parameters have been adjusted Q so that all pa- 
rameters of model F are known. Thus model F is ready 
for application without any further adjustable parame- 
ters. 

In this paper we present a renormalization-group (RG) 
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calculation of the superfluid/normal-fluid interface based 
on model F. The calculation is technically very diffi- 
cult and challenging for two reasons. First the Green 
functions and Feynman diagrams must be evaluated in 
a spatially inhomogeneous system. Secondly, the renor- 
malization factors depend on space and time coordinates 
via the RG flow parameter so that the partial deriva- 
tives with respect to space and time must be replaced by 
appropriate covariant derivatives. 

The first challenge was overcome step by step in several 
previous papers. On the normal-fluid side of the inter- 
face the Green function was calculated [ll| for a zero 
order parameter (ipiz)) — and a linear temperature 
parameter ro(z) ~ + h^z. The local thermal conduc- 
tivity XxiT, Q) and the related temperature profile T{z) 
was calculated. On the superfluid side of the interface 
the Green function and related thermodynamic quanti- 
ties were calculated [l2j for a plane- wave order parameter 
(i^iz)) = jye'*^^ and a constant temperature parameter 
tq. Here a critical superfluid current was found which 
implies a depression of the superfluid transition temper- 
ature T\[Q) < T\ hy & nonzero heat current Q. 

Later the normal-fluid approach [ll| was extended be- 
yond the interface into the superfluid region The 
calculation was made self consistent by a lowest-order 
1 /n expansion which is equivalent to the Hartree approx- 
imation of quantum many-particle physics. In this way 
the superfluid region could be reached where in the whole 
system the average order parameter (tp{z)) = is zero 
due to phase fluctuations related to the motion of vortices 
where however the condensate density ng = {\ip{z)\'^) and 
the superfluid current Jg = (Im['0*V'0]) are macroscopi- 
cally large. The RG theory was applied locally using a lo- 
cal flow parameter which depends on the altitude coordi- 
nate z. The specific heat Cq and the thermal conductiv- 
ity At were calculated for the whole superfluid/normal- 
fluid interface where the effects of the gravity acceleration 
g and the heat current Q were included. The tempera- 
ture profile T{z) was obtained by integrating the heat- 
transport equation Q = —Xt^T. 

In the superfluid region a nonzero temperature gradi- 
ent VT was found which is due to a nonzero thermal 
resistivity induced by the motion of vortices and quan- 
tum turbulence. The theory was especially successful to 
describe the so called self-organized critical state, which 
was predicted by Onuki Q and which was discovered in 
the experiment by Moeur et al. [l5|. In this state the 
temperature gradient VT is equal to the gravity induced 
gradient VTa of pTTj) . i.e. VT = VTa, so that the system 
is homogeneous over a large area in space. 

However, for the superfluid /normal-fluid interface the 
self-consistent approach [l^, works only for large heat 
currents Q }^ Qo — 70nW/cm^ where the heat current 
Q dominates over the effects of gravity g. For smaller 
heat currents this approach does not yield a result. The 
existence and motion of vortices is essential for phase fluc- 
tuations in order to have a zero average order parameter 

mz)) = o. 



For small heat currents Q Qo = 70nW/cm^ vor- 
tices are not present so that the average order parameter 
(-0(2)) is nonzero. In this case, a local calculation is not 
possible. Instead, the full model-F equations must be 
solved as partial differential equations. Here the second 
challenge arises if the RG theory is involved. The RG 
flow parameter is determined locally by a constraint con- 
dition so that it will depend on space and time. This fact 
requires the definition of covariant differential operators. 
A first step for this kind of theory was made by the au- 
thor and Nikodem [l6j . The interface was investigated in 
thermal equilibrium where only the gravity acceleration 
g is present but no heat current. The covariant deriva- 
tives were defined for the renormalized order parameter 
and for the renormalized temperature parameter. The 
renormalized Ginzburg-Landau equation was solved nu- 
merically as a boundary value problem by the multiple- 
shooting algorithm. Results for the order-parameter pro- 
file {il>{z)), the correlation length and the heat capacity 
C were olDtained. However, the calculations [l^ were not 
finished and not published. 

The present paper is devoted to continue, extend, and 
publish our recent calculations [l^. We develop a local 
and time-dependent RG theory for small heat currents 
Q ^ Qo — 70 nW/cm^ in order to fill the gap which our 
previous theory [l^, [111 has left. We solve the partial 
differential equations of model F together with a local 
constraint condition for the RG flow parameter. We cal- 
culate the average order parameter profile {ipiz)) and the 
temperature profile T(z) for the superfluid/normal-fluid 
interface. Furthermore, we calculate the related thermo- 
dynamic and transport quantities, i.e. the specific heat 
Cq and the thermal conductivity At or thermal resistiv- 
ity pt = 1/At- The calculations are not restricted to a 
stationary state of a constant heat current Q. More gen- 
erally, we solve the model-T equations as time-dependent 
partial differential equations, so that time-dependent and 
transient effects can be handled like the propagation of 
second sound. 

The paper is organized as follows. In SecHIlwe briefly 
describe model T, the reduced hydrodynamic model for 
the critical dynamics of liquid "^He near the superfluid 
transition. Furthermore, we explain the approximations 
that we use. In Sec. IIIII we develop our method in order 
to solve the model-T equations together with local con- 
straint conditions for the local RG theory. In Sec. llVl we 
present our numerical results for the superfluid/normal- 
fluid interface in gravity where small heat currents are 
flowing upward or downward. We compare our results 
with the experiment of Chatto et al. p3| and find good 
agreement for the local thermal resistivity. In Sec. fVl we 
compare our small-heat-current results with the large- 
heat-current results of our previous approach 
We discuss the stability of the solutions of our present 
and our previous approach. Finally, in Sec. IVII we com- 
pare our present and our previous approach with other 
theories and recent experiments. We discuss and con- 
clude to which extent our theory can describe mutual 
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friction effects for larger heat currents due to the motion 
of vortices and quantum turbulence. 

II. MODEL AND APPROXIMATION 

The local thermodynamic properties of liquid "^He are 
described by the three standard hydrodynamic variables, 
the mass density p{r,t), the mass-current density j(r, t), 
and the entropy density (j{r,t). Since ^He becomes su- 
perfluid below the critical temperature T\ 2 K , there 
exists an additional fourth hydrodynamic variable, the 
macroscopic wave function ip{r,t), which is the order pa- 
rameter of the superfluid phase transition. The full hy- 
drodynamic equations for superfluid ^He described by 
all these four variables have been derived long ago by 
Pitaevski 

For the critical dynamics near Tx the mass density p 
and the mass-current density j are irrelevant variables, 
because the related hydrodynamic modes, first sound and 
viscosity effects, are fast, p and j can be eliminated or in- 
tegrated out, so that the remaining relevant variables for 
the critical slow modes near the transition (second sound 
and order parameter relaxation) are the order parameter 

and the entropy density a. For these two relevant vari- 
ables, the hydrodynamic equations are given by model F 
and read 

7^ = -2ro7— + ^50^^+0^ . 2.1 
ot oyj* dm 

For convenience and historical reasons, the dimension- 
less entropy density is denoted by m = a/k^. In the 
equations 



H 



d'^r [iro|V^P + i|VVP + Mo|^|^ 

+\Xo^m^ +lam\if -h^m] (2.3) 



is the free energy functional divided by ksT . The Gaus- 
sian stochastic forces O^p and 8m incorporate the fluc- 
tuations. They are defined by the averages {9^) = 0, 
{Om) = 0, and by the correlations 

(e^(r,0^;(r',t')) = 4ro(5(r-r')<5(t-t') , (2.4) 
(0,„(r,t)0,„(r',^')) = -2AoV2<5(r-r')<5(t-<') • (2.5) 

The dimension of the space d is assumed to be arbitrary 
and continuous in the general calculations. However, 
eventually we set d = 3 when evaluating explicit results 
for liquid ^He in a three-dimensional cell. For the cal- 
culations in the critical regime the model-_F equations 
(|2.ip - (|2.3p are treated by field-theoretic means, i.e. per- 
turbation series expansion with respect to Feynman dia- 
grams, renormalization, and the renormalization group. 
For example, the heat capacity and the thermal conduc- 
tivity were evaluated up to two- loop order Q. 



In this paper, we use an approximation following our 
previous work [isl . [l^ . In many-particle physics this ap- 
proximation is known as the Hartree approximation (see 
e.g.Ref.Ei. It is a self-consistent approximation includ- 
ing only a single one-loop diagram which is the tadpole 
diagram. Alternatively, the approximation is obtained 
by the l/n expansion in leading order, where n is the 
number of complex fields in a generalized model with a 
generalized order parameter = {tpi, ...,?/)„). 

For the model-F equations (|2.ip and (|2.2p the approxi- 
mation is obtained by taking the nonequilibrium average 
(. • • ) for all terms and by performing appropriate factor- 
izations of the averages of products of the fluctuating hy- 
drodynamic variables ip, tp* , and m. The factorizations 
are justified by inspection of the Feynman diagrams of 
the Hartree approximation which are shown in Fig. 2 of 
Ref. 0. We factorize the nonlinear terms according to 

2(1^) ~ [T„-V» + 4>io(W.p)+27„(m)](5.) , (2.6) 



5m I 
where 



bm I ' 



(2.7) 



(^)=Xo'M+7o(|^P)-/io . (2.. 



Without an approximation we obtain 



2(Im(V/^)) = V(Im[rVV^]) 



(2.9) 



Consequently, from (|2.ip and (|2.2p we obtain the approx- 
imate equations 



aw 

dt 



d{m) 
dt 



-To [ro - + Auon^ + 270 (m)] {ip) 

+m [xo^{m) + 7o"s - ho] (ip) , (2.10) 

^0 [Xo ^ ('^) + 70^8 - ho] + go VJ^ , (2.11) 



where we define the condensate density Ug and the su- 
perfluid current density Jg by 



= (Im[^*VV^]) 



(2.12) 
(2.13) 



respectively. 

Next, for convenience and simplification of the equa- 
tions we define the temperature parameters 



\ dm I 

= 2X070 [Xo ^("^> + 7ons - ho] , (2.14) 
To = To -I- 2xo7o/io + Aro 

= 7-0 4- 2xo7o [Xo ^("^) + To'^s] , (2.15) 



and the modified temperature parameter 
ri = ro + 4uons , 



(2.16) 
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where uq = uq — ^Xo7o is a combined couphng constant 
following Ref. 0- Thus, the model-F equations can be 
written in the simple form 

M ^ _ro[ri-V2]W+*-^Aro(V), (2.17) 
ot 2xo7o 



d{m) 
dt 



= -Vq 



(2.18) 



The last equation is the heat transport equation where 
(jTi) = ((t) /k-Q is the dimensionless entropy density and 



Ao 



■VAro - .90 Js 



(2.19) 



is the dimensionless entropy current density. The latter 
is related to the heat current Q in standard physical units 
by q = Q/fceT « Q/fceTA. 

The order-parameter equation (|2.17p can be written 
in the form ^("0) = where the operator L is defined 
in (3.13) of our previous paper [ij and related to the 
off-diagonal matrix elements of the inverse Green func- 
tion. This observation shows that the factorizations of 
the present approach are equivalent to the self-consistent 
approximation in our previous paper. We note that the 
factorization is applied only in the order parameter equa- 
tion (j2.17p . The heat transport equation (j2.18p is derived 
without any factorization or approximation. 

The parameters Arp and tq are related to the local 
space and time dependent temperature T = r(r, t) and 
to the critical temperature T\ = Tx{z) of p.ip according 
to 



Aro = 2xo7o = 2xo7o 



Tx 



ro - roc = 2xo7o ■ 



Tx 



(2.20) 
(2.21) 



where Tq is a constant reference temperature. These 
equations have been derived in our previous paper (l^ . 
The critical value of rg is roc = in one-loop approxi- 
mation and hence also in our self-consistent approxi- 
mation. The factor Tx in the denominators is easily ex- 
plained. Since H is the free energy divided by k^Tx and 
since m is the entropy density divided by fee , we find that 
the functional derivative SH/Sm is a temperature divided 
by Tx- We note that the critical temperature Tx = Tx{z) 
defined in (|l.ip depends on the altitude z. Since the gra- 
dient is very small, the z dependence is very weak. Thus, 
in the denominator we may approximately use a constant 
average value which may be the critical temperature at 
the interface z = zq, i.e. Tx = Tx{z) k, Tx{zq). 

Until now, the condensate density Ug and the super- 
fluid current density Jg defined in (j2.12p and (|2.13|) are 
unknown. Since they are defined by an average of two 
fields ip and ip* they are related to the equal-time Green 
function 



G(r,i;r',t) 



(V(r,t))(0*(r',t)) 



(2.22) 



This Green function was evaluated in the Appendix of 
Ref. Q. However, while in our previous paper the aver- 
age order parameter (tp) was zero, in the present paper it 
is nonzero. Hence, we must split the Green function into 
two contributions, a mean-field term and a fluctuating 
term where Sip ^ ip — (ip) is the fluctuating fleld. While 
the mean-field term is expressed in terms of the average 
order parameter {ip) , the fluctuating term is given by the 
result of our previous paper. Consequently, the conden- 
sate density rig and the superfluid current density Jg are 
split into two contributions, too. From Eq. (3.24) and 
(3.25) of Ref. M we obtain 



Us 
J. 



, .90 VAro 1 



l-e/2 



2T'^ 2x070 e 



A„ 



(2.23) 



(2.24) 



Here it is e = 4 — d where d is the dimension of the 
space. Furthermore, Ad = Sdr{l ~ e/2)r(l + e/2) is a 
geometrical factor which is related via Sd — fld/i^ir)'^ 
to the surface of the d dimensional unit sphere ild = 
2 7r''/Vr(d/2). The function $a(X) is defined by the 
divergent series 



r(a- 



N=0 



r(a) 



3A^) X'^ 



(2.25) 



where the argument X is related to the square of the 
gradients of the parameters ri and Aro according to 



X 



12 rf 



(Vri)^ + 2^|(-^VAro).Vri 
To V4xo7oU / 

2n 



90 



4xo7or, 



7 VAro 



(2.26) 



The Green function (|2.22p was evaluated locally for a 
spatial inhomogeneous system where the temperature pa- 
rameters ri and Aro depend on the space coordinate r. 
Gradient terms Vri and VAro are included but curva- 
ture terms and higher derivatives are omitted. This fact 
is clearly seen in the function (|2.25|) and its argument 

Now, all quantities are determined. The approximate 
model-F equations (|2.17p - (|2.18p together with the en- 
tropy current density (j2.19p , the temperature parameters 
(11201), (P^ . (PTTB|) . and the quantities ((^T^ - ip:^ are 
closed equations, which in principle can be solved numer- 
ically. We insert the condensate density (j2.23p into the 
equation for the modified temperature parameter (j2.16p . 
After reordering the terms we obtain 

ri{l + 8uoJ^d$-i+e/2(^)rr^'} =7-o+4uo|(^)P • 

(2.27) 

The left-hand side shows clearly that this is an implicit 
equation for the parameter ri. Furthermore, we insert 
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the superfluid current (|2.24p into the formula for the en- 
tropy current (j2.19p . After reordering the terms we ob- 
tain 



-{ 

o7o 1^ 



1 



9l 



1 



2x070 ^' ■ 2Aor[,£ 
xVAro-.goIm[(V'*)V(7/')] 



(2.28) 



Eqs. (|2.27p and p.28|) of the present paper should be 
compared with Eqs. (3.32) and (3.35) of our previous pa- 
per [iJI , respectively. New contributions are those terms 
on the right hand sides which involve the average order 
parameter {iji). The last term in (j2.28|) may be inter- 
preted as the mean-field contribution of the superfluid 
current. The fluctuating term of the superfluid current 
(|2.24p is proportional to the temperature-parameter gra- 
dient VArQ. For this reason, the fluctuating term is inte- 
grated into the first term of (|2.28p and hence contributes 
to the normal-fluid term. Similarly, in Eq. (|2.27p the 
mean-field contribution of the condensate density is put 
on the right-hand side while the fluctuating contribution 
is put on the left-hand side of the equation. 



III. LOCAL RENORMALIZATION-GROUP 
THEORY FOR PARTIAL DIFFERENTIAL 
EQUATIONS 

The liquid ^He is considered in the critical regime for 
temperatures T close to the superfluid transition at T\. 
In order to treat the critical fluctuations correctly, we 
must renormalize the equations of the previous section 
and apply the renormalization-group (RG) theory. Since 
we consider local physical quantities which are functions 
of space and time, the RG flow parameter will be local 
and depend on space and time. The derivatives with 
respect to space and time in the model- _F equations are 
in conflict with a local RG flow parameter because they 
do not commute with this parameter. For this reason, 
the development of the local RG theory for the model- 
F equations which are partial differential equations is a 
very challenging task. 



A. Renormalization 

We start with the renormalization of the average order 
parameter (V'), the temperature parameters Arg, rp, and 
the coupling constant wq- Following Ref. we have 



Aro 



Ar , 
ZuZ^^in" /Ad) u 



(3.1) 
(3.2) 
(3.3) 
(3.4) 



In these and the following renormalization equations we 
use the convention that the bare quantities are always 
on the left-hand side while renormalized quantities are 



always on the right hand side. The Z factors are the 
renormalization factors. In the Hartree approximation, 
which we use in the present paper and in our previous 
paper [l3|i these Z factors are 



Z^ = 1 



Zr = Zu = 1/[1 - 8u/e] , 



(3.5) 



where it is tqc — 0. The modifled temperature parameter 
ri is not renormalized. We apply the renormalizations to 
Eq. (j2.27|) . multiply both sides with the inverse factor 
Z^^, and reorder the terms. Without any further ap- 
proximation we obtain 



n{l + - 



4u 



Ad 



-l+e/2 



<-'(?)"■"-]} 



/2 



(3.6) 



The average entropy density (m) , the entropy current 
density q, and the remaining model-_F parameters are 
renormalized by Q 

(m) = (xo^m)'/' (mrcn) , (3.7) 

q - (xo^m)'/'qrc„ , (3.8) 

X070 = {xoZm)'^^Zr{fiyAdy/^j , (3.9) 

,90 = {xoZ,ny/^i^iyAdy/^g , (3.10) 

Ao/xo = Z^^ A , (3.11) 

To = Z-^r . (3.12) 

The dimensionless renormalized parameters arc deflned 
by the ratios 



w = r/A , 

F = g/X , 



(3.13) 
(3.14) 
(3.15) 



We note that F = F' -I- iF" and w = w' + iw" are complex 
parameters. The Z factors, which we need explicitly in 
our calculation, arc given in Hartree approximation [l3 | 

by 



Z„,Za = 1/[1 - f/2e] , 



Zr = 1 



(3.16) 



The factor xoZm will cancel out in all our equations. 
Hence this latter factor is not needed explicitly. Applying 
the renormalizations to Eq. (|2.28p we obtain the renor- 
malized heat current 



27 V^' 



X$./2(X) 

1/2 



f_ 

2e 

/2 



2/ 
VAr 



-g (^) Im[(^.;,jV(VA.cn)] . (3.17) 

Again no further approximation is made when reordering 
the terms. In order to evaluate the function $ct(Ar) we 
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need the argument X expressed in terms of the dimen- 
sionless renormahzed parameters. From (|2.26p we obtain 



X = 



12 rl 



(Vri)2+2^f-^VAr) • Vn 
w' \i'yw' I 



w' V47W' 



/ F \2 



(3.18) 



The renormahzation of the model- equations is straight 
forward. From (|2.17l) and (|2.18|) we obtain 

= -r[ri-V2](V'rc„)+*|^Ar(^,,„) ,(3.19) 



dt 



-Vq,.c 



(3.20) 



We furthermore need a relation between the entropy den- 
sity (TOrcn) and the temperature parameters r or Ar in 
renormalizcd form. Wc solve Eq. (|2.15p with respect to 
the entropy density (m), eliminate the condensate den- 
sity Us by (|2.16p . and then perform the renormahzation. 
As a result we obtain 



Ad 



1/2 



27 I 



T 
2u 



1- 



21) 



We have separated the constant value mc,rcn which is the 
entropy at the critical point with temperature T — T\, 
zero heat current Q and zero gravity 5 = 0. We need 
not know this constant value explicitly. Another useful 
quantity is the derivative of (mion) with respect to the 
temperature parameter r. It is related to the renormal- 
izcd specific heat according to 



Ci-c 



27 



\ 1/2 d{m,cn) 



Ad 



dr 



1 



2_ 

2u 



2 



1 



dri 
dr 



(3.22) 



In this way, the time derivative of the renormahzed en- 
tropy density can be expressed in terms of a time deriva- 
tive of a temperature parameter. We find 



dt 



(Adyi^ ar ^ Orcri 

V ^= / dt ~ 27 V / 



1/2 ar a 



27 



'Ad 



1/2 aAr 



dt ■ 
(3.23) 

Since the critical temperature T\{z) does not depend on 
the time, the two temperature parameters r and Ar dif- 
fer by a time-independent value. For this reason, the 
time derivatives of r and Ar are equal. In the present 
paper we prefer the latter time derivative. In this way, 
we reformulate the second model-F equation (|3.20p as 



ac„/AdXi/2 9Ar 



(3.24) 



In the renormahzed specific heat (|3.22l) the remaining 
derivative dri/dr may be obtained as the proportionality 
factor of the gradients Vri and Vr according to 



Vri 



dri 
dr 



Vr . 



(3.25) 



In order to find a relation between the two gradients we 
apply the nabla operator V to Eq. (j3.6p . Thus, we find 

,Vn,{i + f[(i-i)*,,w(IL)--:]} 



Vr + 4M^ V|(i^rc„)|' 
Ad 



(3.26) 



In this result the derivative has increased the index a of 
the function ^a{X) by one. Furthermore, the function 
is multiplied by a factor (1 — e/2). These facts are well 
known from the calculations in our previous paper [l3 |. 
By comparing Eqs. (|3.25p and p.26p we extract dri /dr. 
Since we consider the space dependence only in one di- 
mension z which is the altitude, we obtain a unique re- 
sult. We conclude that in this subsection we have derived 
all equations in renormahzed form which are needed for 
a numerical calculation to solve the model-i^ equations 
as partial differential equations with respect to space and 
time. 



B. Dimensionless renormalized quantities 



In the renormahzation equations (j3.ip - p.4p and p.7p - 
p.l2p the new arbitrary parameter /i occurs which has 
the unit of an inverse length scale. Consequently, this 
parameter may be used to fix the length scale. On the 
other hand in the renormalized model-F equations p.l9p 
and p.24p together with p.l7p the dynamic parameters 
r = r' + zF", A, and g all have the unit of a diffusion 
constant, i.e. length square divided by time. Hence, these 
parameters multiplied by /j,^ may be used to fix the time 
scale. The dimensionless ratios p.l3l) - p.l5p imply that 
only one of these parameters is needed. Thus, we will use 
gfi^ to fix the time scale. 

We rewrite the renormalized model-i^ equations and 
the related renormalized variables and parameters in a 
dimensionless form using fx and g^j? for the scales. Fol- 
lowing our previous paper (l3 | we define the dimension- 
less temperature parameters 



Ap 
P 
Pi 



= Ar/^l^ =T-'{T-To)/Tx 
= r/^i^ = r-i (T - Tx)/Tx , 

= n/p' ■ 



(3.27) 
(3.28) 
(3.29) 



The last equality sign in p.27p and p.28p is obtained by 
renormalizing the bare equations (|2.20p and (|2.2ip . The 
renormahzation factors are combined into the dimension- 
less parameter 



1/2 1 

2^ 



(3.30) 



which may be viewed as a renormalization-group (RG) 
fiow parameter. A change of the length scale by replacing 
fi ^ lii causes a change of r. While £ is the conventional 
RG flow parameter related to the length scale, t = t{£) 
is a RG flow parameter related to the temperature scale. 
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In the hterature on the dynamic RG theory for Uquid ^He 
both flow parameters have been used 0, H, [H, Tl3 . [20j . 
For the dimensionless coupling parameters of model F 
the notations u(£), "f{£), etc. and u[t], 7[t], etc. have 
been used. In the present paper we will use r as the RG 
flow parameter. 

We define the dimensionless renormalized order param- 
eter Y and the dimensionless renormalized heat current 
q by 

(3.31) 
(3.32) 

■ ^ 9 fJ-"^ .90 M" " 
respectively. For convenience of the notation, following 
Ref. [13 we define the dimensionless amplitudes 

A = £-M*-l+e/2(^)pr^'-l] , 

Ai = e-M(l-£/2)<i>e/2(^)pr^'-l] 

Consequently, the renormalized heat current p.l7p can 
be rewritten in the dimensionless simple form 



C. Evaluation of the perturbative amplitudes 



1 



1 



,d-l 



.d~l 



(3.33) 
(3.34) 



|l + ^Ai| VA/9-Im[y*Vy] . (3.35) 



2-fF 



The overall factor l//i is needed to keep the nabla op- 
erators dimensionless. We note that Ad is a geometrical 
factor related to surface of the d dimensional unit sphere 
0. It should not be confused with the amplitudes A 
and Ai. In an analogous way Eqs. p.6p and (j3.26p for 
the modified temperature parameter pi and its derivative 
Vpi can be written in a dimensionless form. We obtain 

p+{Au/Ad)\Y\^ , 
\/p+{Au/Ad)V\Y\^ 



pi{l + 8uA} 
{Vpi){l + 8uA,} 



(3.36) 
(3.37) 



where the second equation should be multiplied by an 
overall factor 1/p to make the nabla operators dimen- 
sionless. Finally, the renormalized specific heat Crcn de- 
fined in (|3.22D and the parameter X defined in p.lSp 
are already dimensionless, so that we can keep them un- 
changed. We must only insert the dimensionless temper- 
ature parameters p.27p - p.29p and use the dimensionless 
nabla operator /i^^V. 

Now, all variables and parameters are expressed in a 
dimensionless form. Thus, we are ready to rewrite the 
renormalized model- equations in dimensionless forms. 
From Eqs. (|3T9)) and ((3^ we obtain 



1 dY _ Wr 



Crcn Ad dAp 

2j gpL^ dt 



+— Apr 

27 



(3.38) 
(3.39) 



where w = 



is a complex parameter. In these 



equations we clearly see that p defines the length scale 
and gp^ defines the time scale. We may interpret p^^V 
as a dimensionless nabla operator and {gp^)~^d/dt as a 
dimensionless time derivative. 



The amplitudes A and Ai , defined in (|3.33p and p.34p , 
respectively, represent the contributions of the perturba- 
tion series expansion which in our case is the Hartree 
term. In order to solve the model- _F equations as partial 
differential equations we must have explicit expression to 
evaluate these amplitudes. The non-trivial contribution 
in the amplitudes is the function $q,(X) together with its 
variable X which are defined in Eqs. (|2.25p and p.lSp . 
This function was first derived in Ref. IllL Unfortunately, 
the function is a divergent infinite series so that it is not 
well defined in this form. However, in thermal equilib- 
rium at zero heat current Q = and zero gravity g = 
this function can be omitted because it is just unity. In 
this case all temperature gradients are zero, so that the 
variable X is zero which implies (f>a{X ~ 0) ~ 1. Hence, 
the amplitudes reduce to 



Ai = e-M(l-e/2)pr/^-l] 



(3.40) 
(3.41) 



Since pi ^ p ^ {T — T\) must be positive, these am- 
plitudes are valid for the normal-fluid equilibrium state 
only and agree with former results 

In the nonequilibrium state the variable X is nonzero. 
In this case the infinite series (j2.25p must be resummed to 
obtain a well-defined expression which can be evaluated 
numerically. Following our previous papers [ill, [3 we 
write 



where C, 



<i>„(x) = [r(a)]-ir-Fa(C) 

-X)-i/3 and 



poo 

T^{Q = / dvv°'-^ exp(-i;3 
^0 



(3.42) 



(3.43) 



The integral is well defined for a > and obtained by 
analytical continuation for a < 0. The new variable C, is 
defined by a third root. Consequently, Q is not unique 
a priori and may be complex. We must specify the root 
which should be taken. For this purpose we define the 
dimensionless parameter 



1 



12, n^^^^)^ 

V47w''^^^J . 



(3.44) 



so that X = —(j/p\. Hence, the new variable can be 
written in the form Q = pi/a^^^ . 

The transition from normal-fluid to superfluid ^He is 
related to a change of sign of pi ~ p ~ {T — T\). Conse- 
quently, also the new variable C changes sign. The non- 
trivial third root to be evaluated is a^^^. For this rea- 
son, we must distinguish two cases which are related to 
the two possible signs of a. This distinction has impor- 
tant physical consequences. There will be two kinds of 
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noncquilibrium supcrfluid phases of liquid ^He which are 
related to the two regimes where either the heat current 
Q or the gravity g is the dominating external influence. 
We discuss these two cases in the following. 



1. Heat current dominated regime: ct > 

The self-organized critical state observed in the experi- 
ment by Moeur et al. implies linear temperature pro- 
files T{z) and Tx{z) as function of the altitude z. The 
temperature difference T{z) — T\{z) = AT is constant 
over a large range of the altitude. Consequently, the re- 
lated gradient parameters Vpi ^ V p ~ V(T — T\) are 
zero. On the other hand the heat current Q causes a 
nonzero constant gradient VAp ^ VT. Thus, in the for- 
mula p.44p only the last term is nonzero which yields a 
positive result for a. Hence, for the self-organized critical 
state the dimensionless parameter a is always constant 
and positive. 

For the inhomogencous noncquilibrium state we may 
conclude that a is also positive whenever the heat current 
Q and hence the related gradient VAp is large compared 
to the effects of gravity. In our previous paper [14| we 
confirm cr > for heat transport in liquid ^Hc on earth 
for heat currents Q > 70nW/cm^. Moreover, for an ex- 
periment in zero gravity conditions in space, a is positive 
for all heat currents. 

Whenever a is positive, the root a^^^ is straight for- 
ward. We just take the real positive root. Consequently, 
the variable C, is real and changes sign at the supcrfluid 
transition. We find C > in the normal fluid regime and 
C < is the supcrfluid regime. The function (|3.42p and 
the integral (j3.43p can be evaluated directly. As a result 
we obtain the amplitudes 



A 



1 

£ 

1 r 

£ 



r(-i + £/2) 
"r(-i + £/2) 



J--l+e/2(C) - 1 , (3.45) 



(3.46) 



which we have derived and used in our previous paper 

We investigate the asymptotic behaviors of the func- 
tion -?q(C) and find 



j-„(C) «r(a)r" 

for C ^ +1 in the normal-fluid regime and 



(3.47) 



F^iO « [i^l'if/'- {-CliT"'"' exp{2(-C/3)-^/n 

(3.48) 

for C <^ — 1 in the supcrfluid regime, respectively |14| . In 
the first asymptotic case p.47p we recover the amplitudes 
p.40p and (|3.4ip of the normal-fluid equilibrium state. 
In the second asymptotic case p.48[) we obtain exponen- 
tially large amplitudes A and Ai for the noncquilibrium 
supcrfluid state. 



The latter asymptotic case has an important physi- 
cal consequence. We consider Eq. (|3.36[) which is a con- 
straint to define pi . In the original form related to (j2.16p 
this equation is rewritten as 



pi= p-^uApi + {Au/Ad)\Y[' 



(3.49) 



The last term is the contribution of the rcnormalizcd 
complex order parameter Y which is nonzero only in 
the supcrfluid state. However, in the supcrfluid regime 
the second term may be a competing term because the 
amplitude A may be exponentially large. Thus, in the 
noncquilibrium system there may be two competing su- 
pcrfluid phases which have different physical properties. 
In Eq. (|2.22p we have split the order-parameter Green 
function into two terms, a mean-field term and a fluctu- 
ating term. The third and the second term in p.49p refer 
to these two terms of the Green function, respectively. 

The complex order parameter Y may be decomposed 
into a modulus rj and a phase according to F = 77 e**^. 
In our previous papers [isl [l^ we argue that in the su- 
pcrfluid regime the modulus 77 and hence the average or- 
der parameter Y is zero because of large fluctuations of 
the phase ip. These large phase fluctuations are related 
to vortices and quantum turbulence. In the present pa- 
per we consider a nonzero average order parameter Y in 
the supcrfluid regime and solve the renormalized model- 
F equations numerically as partial differential equations. 
We find a competition between the mean-field superfluid 
phase, described by the average order parameter Y , and 
the fluctuating superfluid phase, described by the expo- 
nentially large amplitude A. 



2. Gravity dominated regime: a < 

In thermal equilibrium for zero heat currents Q = the 
temperature T is constant. Consequently, the gradient 
VAp ~ VT is zero. On the other hand, gravity on earth 
implies a nonzero gradient of the critical temperature 
VTa- Hence the other gradients Vpi~ Vp ~ V(r - 
Tx) = -VTa are nonzero. In Eq. (j?^ only the first 
term is nonzero which implies a negative dimensionless 
parameter a. A small heat current Q will not change the 
situation. In our numerical calculations we find ct < 
for Q < 20nW/cm2. 

An exception is the self-organized critical state which 
always implies cr > and which exists for arbitrary small 
heat currents Q where the temperature difference T(z) — 
T\{z) = AT is constant. Nevertheless, for small heat 
currents AT is positive so that the system is normal fluid 
and the sign of a is irrelevant. 

For negative a the third root is always complex. We 
flnd CT^/^ = e^"/^(— cr)^/^, so that the variable of the 
function ((Hlig)) is complex, i.e. C = e'^'''^ 1'^ p\l {-ofl'^ . 
For convenience we introduce the new real parameter 
C = P\I{—^Y^^ which is related to the old parameter 
via C = e^^'^/^C- We furthermore deflne the new complex 



9 



function 

GaiC) = 6^*""/=^ J-„(C) = 6^*""/=^ -F«(e±*"/'C) • (3.50) 

which can be decomposed into real and imaginary parts 
according to GaiC) = G'aiC) ^iG'aiC)- The new complex 
function is not uniquely defined because there are two 
complex roots which can be chosen. This fact causes 
two possible signs for the imaginary part. However, we 
choose the so called principal part, which is obtained 
as the average of the two cases so that the imaginary 
part cancels. Thus, we simply omit the imaginary part 
GaiO- In the normal- fluid region C > this assumption 
is plausible because the imaginary part converges to zero 
exponentially for increasing C. As a result, we rewrite the 
amplitudes p.45p and (|3.46p in terms of the new function 
(|330)) as 



A 



{-<y) 



-e/6 



r(-l + e/2) 



r(-i + e/2) 



r'e-i+e/2(c)-i 



G'edO-^ 



(3.51) 
(3.52) 



Once again, we consider the asymptotic behaviors of the 
function G'aiO- We find 



ea(C)«r(«)r" 



(3.53) 



for C +1 in the normal-fluid regime and 

G'aiO « (7r/3)i/2(-C/3)"/2-3/4 

X cos{2(-C/3)3/2 ^ (7r/4)(2a - 1)} (3.54) 

for C <^ — 1 in the superfluid regime. In the first asymp- 
totic case p.53p we recover the amplitudes p.40p and 
p.4ip of the normal-fluid equilibrium state. Since here 
the amplitudes do not depend on the dimensionless pa- 
rameter a at all, in the normal-fluid regime the sign of 
a is irrelevant. In the second asymptotic case (|3.54p the 
function ^q(C) and hence the amplitudes A and Ai os- 
cillate but remain of order unity. 

Again, the latter asymptotic case has an important 
physical consequence. In Eq. p.49|) the second term is al- 
ways small because the amplitude A never becomes large. 
Hence the superfluid phase is unique. It is the mean-field 
superfluid phase where the average order parameter Y is 
nonzero. Vortices due to fluctuations effects and a fluc- 
tuating superfluid phase do not exist for cr < 0. 



D. Renormalization-group theory 
and flow parameter condition 

In the renormalization procedure the parameter /i is 
introduced which fixes the length scale. This parame- 
ter generates a transformation group which is known as 
the renormalization group. Following Ref. it can be 
changed by the substitution ^ fi£, where the dimen- 
sionless parameter £ is called the renormalization-group 



(RG) flow parameter Q. However, for simplicity and 
consistency of the following calculations, in this paper 
we do not use the above substitution. We avoid the use 
of the flow parameter £ and thus change the length scale 
parameter fi directly. We use the alternative dimension- 
less RG flow parameter r which is defined in p.30[) . All 
quantities of the renormalized theory can be expressed 
in terms of this RG flow parameter. The dimensionless 
coupling parameters are u[t], 7[t], w[t] — w'[t] 4-m"[T], 
-F[t], and fir]- This is a notation which was defined in 
Refs.Qandi. 

A differential relation between the flow parameter r 
and the length-scale parameter ^ can be obtained by a 
logarithmic differentiation of Eq. p.30[) , which reads 



dlnr 



1 



d\nZ- 



9 In 7 



- 2 V 9 In /i / 9 In /.t 
Using the definitions of the RG zeta functions 

= ainZ^V^ln/i, 
Cr = dlnZ-'^/d\nfi , 
Cm = 91nZ,7//ain/i , 

using the RG equation for the parameter 7 01 
91n7/ainAi- [-£ + 2C. + Cm]/2 , 



dln^i . (3.55) 



(3.56) 
(3.57) 
(3.58) 



(3.59) 



and using e = 4 — d, 



p.55p can be simplified into 



ff In T = [2 — Cr] din ^J. 



(3.60) 



The zeta function Cr = Cr(u) is explicitly available as a 
function of u = u[t] [21|. Thus Eq. p.60p enables an 
explicit numerical calculation of r as a function of ^ and 
vice versa. 

Since the renormalization procedure implies a reorder- 
ing of the perturbation series, the RG flow parameter t 
should be chosen in an optimum way so that the conver- 
gence behavior of the series is optimized. To do this we 
choose the constraint condition 

3pi -2p + 3{Au/Ad)fY(yY/fif + fAp(yAp/^l)^ = 1 . 

(3.61) 

The modified temperature parameter pi is defined in 
Eq. p.36p which may be viewed as a second constraint 
equation. The first two terms on the left-hand side of 
p.6ip guarantee the standard flow parameter conditions 
of normal-fluid and superfluid ^He in thermal equilibrium 
and zero gravity which have been formulated in Ref. 0. 
The latter two terms are gradient terms which stabilize 
the intermediate region of the superfluid/normal- fluid in- 
terface. The two parameters fy and /ap are dimension- 
less and control the influence of the gradient terms. In 
our calculations we have used /y = 5 and Jap = 1 as an 
optimum choice. 

In thermal equilibrium and zero gravity all quantities 
and parameters are constant in space and time. An ex- 
ception is the renormalized order parameter Y = Y(t) = 
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^^%v(t) together with the constant modulus rj and the 
time-dependent phase ip{t) = —ijjt + tpQ. Since all gradi- 
ent terms are zero, the model-F equations (|3.38p . p.39p 
and the flow-parameter equation p.6ip reduce to 



Pi y = , 

3pi - 2p = 1 . 



(3.62) 
(3.63) 
(3.64) 



The first equation is always satisfied because it defines 
the order-parameter frequency lo in terms of the dimen- 
sionless renormalized temperature difference Ap where 
the time scale is ruled by the parameter combination 

In the normal-fluid state, the second equation p.63p 
implies the zero order parameter y = 0, where pi may 
be nonzero. The flow-parameter equation p.64p together 
with the constraint p.36p and the amplitudes p.40p and 
([?:iT|) imply p = pi = 1, ^ = 0, and Ai = -1/2. These 
results are compatible with the equilibrium theory of Ref. 
0. The resulting flow parameter condition is p = 1 which 
in the notation of Ref. reads r{l)/{^if = 1. Conse- 
quently, from Eq. (|3.28p we obtain the flow parameter 
T ^ {T — T\)/T\ which just is the reduced temperature 
as known from earlier work 

In the superfluid state, Eq. p.63p implies pi = 
where the order parameter Y is nonzero. Consequently, 
Eq. p.64p yields the flow-parameter condition p = 
— 1/2 which is well known from Ref. in the nota- 
tion r(€)/(p£)2 = -1/2. Again, Eq. ((X^ relates the 
flow parameter to the reduced temperature according to 
r = 2{T\ — T)/T\. From the constraint p.36p wc obtain 
the modulus of the order parameter -q = \Y\. Since the 
left-hand side is zero, we obtain 77 = {Adl?>uY^'^ . 

The above investigation of the normal-fluid and super- 
fluid equilibrium states in zero gravity shows, that in our 
numerical calculations for the superfiuid/normal- fluid in- 
terface the dimensionless renormalized temperature vari- 
ables p, pi, and the modulus of the dimensionless renor- 
malized order parameter = \Y\ must approach constant 
asymptotic values on both sides far away from the inter- 
face. In the intermediate region near the interface, the 
variables will interpolate the asymptotic values. The RG 
fiow-parameter condition (|3.6ip guarantees the asymp- 
totic values and yields an appropriate interpolation in 
the intermediate interface region. The gradient terms in 
this condition will stabilize the interpolation. 

The RG flow-parameter condition p.6ip is designed for 
the superfiuid/normal- fluid interface at small heat cur- 
rents where gravity is the dominating external influence 
and where in the superfluid phase the order parameter Y 
is nonzero. In our classification of Sec. IIII CI this super- 
fiuid phase is the mean-field superfluid phase. The other 
case is the fluctuating superfluid phase where the order 
parameter Y is zero and vortices are present. This latter 
case has been investigated in our previous publications 
(13L [l4| where the RG flow-parameter condition is given 
by Eqs. (11) and (4.39) of Refs. [H and [11 respectively. 
This latter flow parameter condition can be compared 



with our present condition p.6ip if p is eliminated by us- 
ing the second constraint (j3.36p and if we use pi = cr^/^C. 
Then, the first and second term of our present condition 
p.6ip are identified with the second and third term in 
Eqs. (11) and (4.39) of Ref. [H and [ll. The gradient 
terms of Eq. p.6ip are replaced by the first term in Eqs. 
(11) and (4.39) of Refs. [l3l and Il4l. which is also a gradi- 
ent term because a is depends on the gradients following 
(13.441). We note that the RG fiow-parameter condition of 
Ref. I13I and [H is designed for the fiuctuating superfluid 
phase where the order parameter Y is zero and vortices 
are present. 



E. Covariant derivatives 



The flow parameter equation p.6ip and the constraint 
condition p.36p are local equations. Consequently, the 
flow parameter r, the renormalization Z factors, and 
the dimensionless coupling parameters are local and de- 
pend on space and time. This fact will affect the space 
and time derivatives in the renormalized equations. We 
must replace the partial differential operators by covari- 
ant derivatives. To do this, we write the renormalization 
equations in a form like Eqs. p.ip - p.3p . so that the bare 
quantities are on the left-hand side and all renormalized 
quantities are on the right hand side. Then we apply the 
differential operator. We start with the renormalization 
of the temperature p.27p which is equivalent to p.2p . 
We apply the nabla operator and obtain 

V[(r-To)/TA] = V[tAp] =T[V + (Vlnr)]Ap 

= rDAp . (3.65) 

The last equality sign defines the covariant derivative. 
We continue with the temperature difference p.28p which 
is equivalent to p.3p and proceed in the same way. As a 
result wc obtain the covariant derivatives 

DAp = [V + (Vlnr)]Ap , (3.66) 
Dp = [V + (V In r)]p . (3.67) 

The covariant derivative of pi is more complicated. Gen- 
eralizing Eq. p.37p we obtain 



(Dpi){l + ^uAi} = Dp + i^u/Ad) D|y|' 



(3.68) 



which can be resolved with respect to Dpi . 

Next we consider the renormalization of the order pa- 
rameter p.ip . We write this equation in terms of the 
dimensionless renormalized order parameter Y by using 
p.3ip . apply the nabla operator, use p.56p . and obtain 



^l/2^(<i-2)/2 
^l/2^(,_2)/2 



1 



2 
Id 

2~ 



(d-2 



2-C, 



^(Vlnr) 



^l/2^(d-2)/2Dy 



y 
y 

(3.69) 
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Again, the last equality sign defines the covariant deriva- 
tive. We define the running critical exponents [2l| 



V = 1/(2 - C.) , 
/3 = zy(d-2 + ?7)/2 



(3.70) 
(3.71) 
(3.72) 



These exponents are called running exponents because 
they depend on the RG flow parameter r via the zeta 
functions and thus carry all the Wegner corrections. In 
the asymptotic limit t ^ they converge to the univer- 
sal critical exponents. Then from Eq. p.69p we obtain 
the covariant derivative of the dimensionless renormal- 
ized order parameter 



Dy = [V + /3(Vlnr)]y 



(3.73) 



We note that /? is the running critical exponent of the 
order parameter. This result makes clear, how the gen- 
eral structure of a covariant derivative of a dimensionless 
renormalized quantity looks like: It is the partial deriva- 
tive of the quantity plus the critical exponent times the 
partial derivative of In t times the quantity. In the model- 
F equations we also need the second covariant derivative 
of the order parameter. It is obtained by applying the 
operator twice, i.e. 



D^y = [V + /3(V In r)]2y 



(3.74) 



Furthermore we consider the renormalization of the 
heat current p.32p . Applying the nabla operator we ob- 
tain 



V[Q/5oA:bTa 



/^MV + (rf-l)(Vln/i)]q 

[V + (rf- l)!^(VlnT)]q 
/-^Dq. (3.75) 



Thus, we find the covariant derivative of the dimension- 
less renormalized heat current 



Dq= [V + (rf- l)i^(Vlnr)]q 



(3.76) 



We identify {d — \)v as the running critical exponent of 
the heat current. The inverse exponent x — l/[(rf — l)t^] 
is known from the depression of the critical temperature 
T\ by a nonzero heat current Q 0, [l^l . We note that 
the covariant derivatives ^^-^^ and ([X75)) - ([5:7I)) 
have been derived already in our previous unpublished 
approach [l6| for the interface in thermal equilibrium at 
zero heat current. 

Above, we have defined the covariant derivatives with 
respect to the space coordinates D. We also need the 
covariant derivatives with respect to time Dt- To obtain 
them we replace the nabla operator by the partial time 
derivative dt = d/dt. Thus, as results we obtain e.g. 

DtY = [a* -f/3(at In r)]y , (3.77) 
AAp = [St-f (9tlnr)]Ap . (3.78) 



Now, we are ready to rewrite the model-F equations 
in terms of covariant derivatives. From Eqs. p.38p and 
(13.391) we obtain 



27T 50 ^ 



27 



Apy 



(3.79) 
(3.80) 



In these equations wc have performed some further sub- 
stitutions which are known from our previous paper |14| , 
i.e. 



M = r 



(3.81) 



Here ^ = ^[r] is the correlation length. Close to criticality 
it has the asymptotic form ^ = The identifica- 

tion ^ — IS correct in our Hartree approximation 
which is a self-consistent one-loop approximation. Cor- 
rections appear in higher orders [2l|. The renormalized 
time-scale parameter g^^ is expressed in terms of the 
bare parameter 50 by using the renormalization equation 
p.lOp where XoZm has been eliminated in favor of t by 
(|3.30p . As a result the renormalized model-F equations 
(|3.79p and (|3.80p are dimensionless equations for the di- 
mensionless quantities. There are two parameters, which 
control the scales of space and time. They are ^0 and go, 
respectively. 

The model-i^ equations (j3.79p and p.80p are supple- 
mented by some further equations including the dimen- 
sionless renormalized entropy current (|3.35p and the two 
constraint conditions (|3.36p and (|3.6ip where all nabla 
operators V are replaced by respective covariant deriva- 
tives D. Thus, we obtain the dimensionless renormalized 
entropy current 



A, 



[l + ^A,}^T>Ap-lm[Y*^BY] , (3.82) 



27F 
the constraints 



3pi-2p + 3(4uMd)/y(eDy)2 
+/Ap(eDAp)2-l = , 
K2 - pi{l + 8uA}-[p+{4u/Ad)\Y\^] = 

and furthermore the dimensionless variables 



(3.83) 
, (3.84) 
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(eDpi)2 + 2- 



F 



;^DAp -(eDpi) 



F 



47W 



(3.85) 



and ( = pxja^l'^ or C, = pi/(— tr)^/'^ which are needed to 
calculate the dimensionless amplitudes (j3.45p - (j3.46p or 
(imD-(l53a. 
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F. Numerical algorithm 



The numerical algorithm for solving the renormalized 
model-F equations p.79p and 
constraints (j3.83p 
ations. First on the left-hand sides of the model-F equa- 
tions the partial time derivatives within the covariant 
derivatives arc replaced by discrete forward differences 



p.80[) together with the 
and p.84p is implemented by two iter- 



dtY 
dtAp 



[Y{r,t + At) - 
[Ap{r,t + At) 



Y{r,t)]/At, (3.86) 
-Ap{r,t)]/At. (3.87) 



Secondly, the constraints are solved by a Newton method. 
The two iterations are performed in parallel, i.e. alterna- 
tively one time step and one Newton step. In this way 
starting with appropriate initial functions at an initial 
time the dimensionless renormalized quantities Y{r,t), 
Ap(r,i) and pi{r,t), \nT{r,t) are obtained as functions 
of space and time. All the other dimensionless renormal- 
ized quantities which are needed on the right and sides 
of the iteration equations can be calculated from the four 
quantities by formulas we have derived above. The co- 
variant derivatives DF, D^F, DAp, Dp, and Dq are 
calculated with discrete nabla and Laplace operators on 
an equidistant grid of the space-coordinates r. The co- 
variant derivatives of further quantities can by related to 
those five by equations like (|3.68|) . 

For the Newton-iteration step we need the derivatives 
of the constraint functions Ki and K2 with respect to pi 
and Inr. We use the derivatives 



dKi/dpi 
8X2/ dpi 
dKi/dlnr 

3X2/ din T 



3 , (3.88) 
1 + 8uAi , (3.89) 
-{-2p + 3(4uMd)/y(^Dy)22(/3 + y) 
+/Ap(eDAp)22(l + j/)} , (3.90) 
-{MEi2{l + v) 

-[p+{Au/Ad)Y^2p]] (3.91) 



together with the amplitude 

El = (pi/6)[(2-e)^-2^i -1] 



(3.92) 



The latter two derivatives are approximations, because 
wc omitted the weak dependence of the dimensionless 
coupling parameters u[t\^ 7[t], etc. on the logarithmic 
RG flow parameter Inr. Nevertheless, our numerical cal- 
culation works. There is no significant influence of this 
approximation. 

Our numerical calculations are performed very close to 
criticality where r < 10~^. Consequently, for the running 
exponents we can use the universal critical exponents as 
a good approximation. We use the experimental value 
V = 0.671 of Lipa et al. [l^, and the theoretical value 
77 = 0.038 of Schloms and Dohm The exponent /3 = 
0.348 is calculated from the scaling relation p.72p . where 
the dimension of space is = 3. Finally, we can use the 
asymptotic formula for the correlation length ^ = '''"'^ 
as a good approximation. 



Our numerical calculations show that the iterations are 
stable for small heat currents Q = |Q| ^ 20nW/cm^ 
where the gravity is the dominating external force and 
the dimensionless parameter a defined in (|3.85p is al- 
ways negative. For lager heat currents the parameter a 
will have a sign change locally in space, which causes nu- 
merical troubles. We can stabilize the calculations up to 
a maximum heat current Qmax = 160nW/cm^ by adding 
a small imaginary constant to the right-hand side of Eq. 
(j3.85p . However, for larger heat currents where the heat 
flow is the major and the gravity is the minor external 
influence the iteration is unstable so that no results can 
be obtained. 



IV. NUMERICAL RESULTS 

Most experiments with liquid ^He close to the super- 
fluid transition are performed at saturated vapor pres- 
sure. The temperature T is varied in the region near T\ 
where the pressure is kept at the value of the liquid- 
gas transition. In this case the critical temperature 
is Tx = 2.172 K. The parameters which specify the 
scales of length and time are ^0 = 1-44 x 10~^cm and 
go = 2.164 X 10"s-\ respectively [13, [H, [H] . The di- 
mensionless renormalized coupling parameters u[t], 7[r], 
w[t\ = w'[t\ -\- iw"[T\, F[t\, and /[r] as functions of the 
RG flow parameter r are taken from Ref. H. 

We perform the numerical calculations for liquid ^He in 
d = 3 dimensions. The system is assumed to be homoge- 
neous in the two horizontal directions x and y. Thus, all 
quantities and functions depend only on the altitude co- 
ordinate z and the time t. The model- equations reduce 
to partial differential equations with the two variables z 
and t. The size of the experimental cells which contain 
the liquid ''He is usually some millimeters in z direction. 
We use a cell length L = 2.0 mm and discretize the z 
coordinate into 500 points. Consequently, the discretiza- 
tion is A2 = 4.0 pm in the altitude coordinate. 

The discretization of the time At in the partial deriva- 
tives (|3.86p and p.87p must be sufficiently small so that 
the iteration converges. On the other hand At should 
be sufficiently large, so that the calculation time on the 
computer is not too long. We find At = 4.0 x 10~® s 
as an optimum choice. Starting the calculations in any 
nonequilibrium state, we first observe space and time de- 
pendent oscillations which are related to second sound. 
These oscillations relax on a time scale of about one sec- 
ond. After a time interval St = 2.0 s the system reaches a 
stationary state with a constant homogeneous heat cur- 
rent Q where all oscillations are disappeared. This means 
we need 5 x 10^ iteration steps on the computer until the 
system converges to the steady state. 

For a heat flow in z direction there must be a heat 
source and a heat sink at the boundaries of the cell 
zi = —L/2 = —1.0 mm and 22 = +L/2 = +1.0 mm, re- 
spectively. Thus, a source and sink term must be added 
to the heat transport equation of model F (|2.2p which is 
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given by 

W{r, t) = 2 [Qi 5{z - zi) - Q2 5{z - Z2)] 



(4.1) 



In dimcnsionless renormalized form the source and sink 
term is 



w 



(r, t) = 2 [gi 5{[z - z,]/i) - q2 S([z - Z2]/0] ■ (4.2) 



This latter term must be added to the dimensionless 
renormahzed model- equation p.80|) on the right-hand 
side. The relation between the dimensionless renormal- 
ized heat currents q and the physical heat currents Q 
is obtained from the renormalization equations 
(Eini), and (13321) • We obtain 



(4.3) 



50 9okB.T\ 



which should be applied to both heat currents in (|4.ip 
and (|4.2p . It is important to note, that the RG flow 
parameter r = T(r,i) and hence the correlation length 
^ = ^[r] depend on space and time. This fact is impor- 
tant for Eqs. (|4T1) and (|42l) . 

We perform the calculations in the following way. First 
the system is stabilized in the thermal equilibrium. Then 
at time t ^ to the external heat source and sink (|4.ip 
or (|4.2p is switched on where we chose equal values 
Qi ^ Q2 = Q- Then after a time interval 6t = 2.0 s all 
oscillations are relaxed and the system reaches a steady 
state. The local heat current Q(r,i) = Qbz will be ho- 
mogenous in space, constant in time, and directed verti- 
cally along the z axis. 

The boundary conditions at zi and Z2 are important 
for the stability of the iterations. There should be no 
boundaries at all. This means we need periodic bound- 
ary conditions. The system can be made periodic in the 
following way. We mirror the cell at one of the bound- 
aries. Then wc obtain a periodic structure of length 2L. 
Furthermore, for the discretization the delta functions in 
Eqs. (j4.1|) and (j4.2p must be replaced by smooth peaks 
of a small width Sz. We choose 5z = 3 Az which is a few 
discretization lengths. From the heat source at 2 = zi the 
heat current Qi will flow away in both directions, where 
on the other hand a heat current Q2 will flow from both 
directions to the heat sink at Z2. This fact explains the 
factor 2 in Eqs. (|i?T|) and 



A. Dimensionless renormalized quantities 

The direct results of the numerical calculation are the 
dimensionless renormalized temperature parameters Ap, 
p, pi, and the dimensionless renormalized order param- 
eter Y as functions of the altitude coordinate z and 
the time t. In Fig. [1] the results are shown for the 
superfluid/normal- fluid interface of liquid ^He at zero 
heat current Q = in thermal equilibrium. The interface 
is induced by the gravitational acceleration g = 9.81 m/s^ 
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FIG. 1: (Color online) The dimensionless renormalized tem- 
perature parameters p (green solid line), pi (blue dashed line), 
and the modulus of the order parameter 77 = |y| (red dash- 
dotted line) for the superfluid/normal-fluid interface. 



on earth. Since in thermal equilibrium the temperature is 
constant we may choose it equal to the reference temper- 
ature so that T = To- Hence Eq. (|3.27|) implies Ap = 0. 
This is a trivial result which is shown by the black dot- 
ted line. The parameter p is related to the tempera- 
ture difference T — T\ by p.28p and shown as green solid 
line. The modified temperature parameter pi is defined 
in p.36p and shown as blue dashed line. Finally, the 
modulus of the dimensionless renormalized order param- 
eter 7] = \Y\ is shown as red dash-dotted line. 

In Fig. [T] we observe three different regions. For low 
altitudes z < — 100/im we find the asymptotic values 
/9 — > 1, pi — !> 1, and 7] = \Y\ — !> 0. Hence, in this re- 
gion the "'He is normal fluid. We recover the related 
flow parameter condition p = 1 of Ref. in the asymp- 
totic limit z — > —00. For high altitudes z > -1-100 /^m 
we find the asymptotic values p — > —1/2, pi 0, 
and ri = \Y\ ^ (Ad/Su)'^/'^ where Ad = l/47r for 
d ~ 3. Hence, in this latter region the order parameter 
is nonzero and the ^He is superfluid. Again, we recover 
the related flow parameter condition p ~ —1/2 of Ref. 
in the asymptotic limit z +00. The third region is 
the interface region —100/im < z < -1-100 pm. Here the 
curves interpolate between the asymptotic values. We 
clearly see that the interface induced by gravity has a 
thickness of about Azj g « 200 p-m. 

Since the system is constant with respect to the hor- 
izontal coordinates x, y, and with respect to the time 
t, the covariant derivatives of the dimensionless renor- 
malized quantities are nonzero only for the altitude co- 
ordinate z. In most cases these covariant derivatives are 
calculated by numerical differentiation using the formu- 
las derived in section UlI El An exception is ^.D^pi which 
is expressed in terms of other covariant derivatives by 
formula p.68[) . The result is shown in Fig. [5] by the blue 
dashed line. Alternatively, we apply Eq. p.67p to the 
modifled temperature parameter pi and calculate the co- 
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FIG. 2: (Color online) The blue lines show the dimensionless 
covariant derivative ^D^pi calculated in two ways: by nu- 
merical differentiation (blue solid line) and by formula (|3.68p 
(blue dashed line). Furthermore the parameter a defined in 
(|3.85p and multiplied by a factor 12 is shown as orange dash- 
dotted line. Finally, the argument of the function 1)3. 50|) ( is 
shown as magenta double-dash-dotted line. 



variant derivative directly by numerical differentiation. 
This latter procedure is not correct in the interface re- 
gion where the renormalization factors depend on the 
altitude coordinate because pi is not renormalized as p. 
Nevertheless, in Fig. [2] the result is shown by the blue 
solid line. Surprisingly, the two blue lines, the solid one 
and the dashed one, are close to each other. Hence, Eq. 
(|3.67|) is not that bad for calculating the covariant deriva- 
tive ^D^pi. 

The blue lines in Fig. [2] represent the covariant deriva- 
tive of the blue dashed line in Fig. [T] However, the latter 
line represents pi = pi{z) and shows a negative min- 
imum value pi.min < at the position ^min ~ 60 /.im. 
Consequently for S^D^pi we expect a zero at this position 
related to a sign change. In Fig. [3] the solid blue line does 
show this zero and sign change but the dashed blue line 
does not. In this way, the apparently incorrect formula 
(j3.67p for pi appears to be more realistic than the generic 
formula 

The existence of the sign change is supported by the 
following argument. In a small z interval close to the 
interface we may modify the renormalization-group the- 
ory by choosing a constant flow parameter r. In this case 
the covariant derivatives reduce to the partial derivatives 
so that Eqs. p.67p and p.68p would yield identical re- 
sults for ^DzPi ~ £,dzPi and the two blue lines in Fig. 
[2] would collapse to a single line. As a result the sign 
change would be found at Zmin if we evaluate the partial 
derivative explicitly by differentiation of the blue dashed 
line in Fig. \T\ 

However, the sign change of the solid blue line in Fig. 

would have a dramatic consequence for the numerical 
procedure when calculating a and the amplitudes A and 
Ai . In thermal equilibrium we have Ap = so that Eq. 



p.85p reduces to ct = —{£_DzPi)'^/12. Consequently, a 
will be negative everywhere except at a point close to 
Zniin. At this point we have cr = so that the formulas 
for the amplitudes A and Ai reduce to (|3.40p and (|3.4ip . 
respectively. However, close to the minimum position 
^min the modified temperature parameter pi is negative 
which implies an imaginary result for p^ "^^ in Eqs. p.40p 
and (|3.4ip where e = 4 — d = 1 for = 3. Hence, the 
amplitudes A and Ai are not well defined if the solid line 
in Fig. [2] and the formula p.67p is used. 

The problem arises due to the fact that we evaluate the 
Green function (|2.22|) . the condensate density ng, the su- 
perfluid current Jg, and hence the amplitudes A and Ai 
in an approximation where only the covariant gradients 
of the temperature parameters pi and Ap are taken into 
account. If we could do the calculation for the full space 
dependence all these quantities would be well defined. In 
the unpublished work [iGj we extended the calculation 
by including also the curvatures of the temperature pa- 
rameters. While the problem at z,nin was abolished, the 
calculation was much more complicated and restricted 
to the thermal equilibrium at zero heat current. More- 
over, other mathematical difficulties appeared. Hence, 
this more sophisticated calculation could not be realized 
in practice for our purpose. 

However, our fortune is a small inaccuracy of the 
approximation in our numerical calculation in practice 
which implies that the blue dashed line in Fig. [5] is com- 
pletely negative and does not show a zero and a sign 
change for the covariant derivative £.DzPi. The related 
parameter a ~ — (^Dz/9i)^/12 is shown in Fig. [5] as or- 
ange dash-dotted line where it has been enhanced by a 
factor 12. Clearly, this curve is negative and never zero 
for all altitudes z in the interface region. For this reason 
we can apply our formulas (|3.5ip and (j3.52p for the am- 
plitudes A and Ai without a problem if we use the generic 
formula p.68p for the dimensionless gradient S^D^pi- We 
obtain smooth and stable results which are within the 
accuracy of our approximation. 

In order to evaluate the amplitudes A and Ai we need 
the function (|3.5Gp and its argument C, = cr)^/'^. 
Consequently, from the dashed blue line in Fig. [T] and 
the orange dash-dotted line in Fig. [5] we obtain the di- 
mensionless variable ^ as a function of the altitude coor- 
dinate z which is shown in Fig. [2] by the magenta double- 
dash-dotted line. In the normal-fluid region for z < 
the variable ^ increases quickly for decreasing altitude z. 
Consequently, in this case the asymptotic formula p.53p 
can be used so that the amplitudes A and Ai reduce to 
the simple formulas (|3.40p and (|3.4ip of the normal-fluid 
equilibrium state. In the superfluid region near the in- 
terface the variable C is negative. However, it is bounded 
from below by the value —1. Consequently, the asymp- 
totic formula (|3.54p is not needed. This means that the 
variable C never comes in the large negative region where 
the function (|3.50p oscillates and possesses a significant 
imaginary part. This observation is very important for 
the consistency of our theory because the oscillations 
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would be unphysical and the imaginary part would be 
related to an instability. 



B. Temperature profiles 

Until now, the calculations are restricted to the ther- 
mal equilibrium at zero heat current Q = 0. Here the 
phase of the order parameter Y = rje^'^ is constant, so 
that wc can chose (p ~ 0. We have extended our nu- 
merical calculations to small nonzero heat currents Q in 
the interval -TOnW/cm^ <Q< +160nW/cm2. In this 
latter case the phase of the order parameter ip = ip(z,t) 
will be a nontrivial function of the altitude coordinate z 
and the time t. A positive heat current Q > means a 
heat flow Q ~ Q in the z direction which means that 
the heat current flows upward from bottom to top. The 
original experiment by Duncan et al. and succeeding 
experiments investigating the superfluid/normal-fluid in- 
terface induced by a heat current Q were performed in 
this configuration. On the other hand a negative heat 
current Q < means a downward heat flow from top to 
bottom. This latter configuration was investigated much 
later in the experiment by Moeur et al. [Tsj . 

In the nonequilibrium system with a nonzero heat flow 
the dimensionless rcnormalized temperature parameter 
Ap will be nonzero. Once the local space and time depen- 
dent RG flow parameter r = T{z,t) is known, the space 
and time dependent temperature profile T = T{z,t) is 
calculated from Ap by Eq. p.27p . Furthermore, the lo- 
cal space and time dependent heat current Q = Q{z,t) 
is calculated from the dimensionless renormalized heat 
current q by Eq. (|4.3p . After a time difference of about 
6t ~ 2s the system will relax in a stationary state where 
all quantities are constant in time. If in Eq. (|4.ip the 
source and sink parameters Qi = Q2 = Q arc chosen, a 
vertical heat current Q will be found in the whole system 
which is constant in the space variable z. Consequently, 
the different temperature profiles we obtain in our nu- 
merical calculations can be labeled by this constant heat 
current. 

Our numerical results are shown in Fig. El The tem- 
perature profile T{z) is shown by the colored solid lines 
for several values of the heat current Q which are speci- 
fied in the caption of the figure. On the other hand, the 
superfluid transition temperature Tx{z) as a function of 
the altitude coordinate z is shown by the straight black 
dashed line. The slope of this latter line is the effect of 
the gravity on earth. 

The altitude zq at which the temperature profiles T{z) 
and Tx{z) intersect each other so that T{zq) ~ T\{zq) 
may be viewed as a reference altitude to specify the posi- 
tion of the superfluid/normal-fluid interface. We realize 
that the system is translation invariant in the sense that 
we can move the curves parallel along the straight dashed 
line. Thus, for convenience and simplicity we select a co- 
ordinate system so that all curves intersect at the same 
altitude z = zq = 0. This choice is no physical restriction 




400 



FIG. 3: (Color online) The temperature profiles T{z) of the 
superfluid/normal-fluid interface of liquid '*He in gravity are 
shown for several heat currents Q as colored solid lines. The 
solid lines on the left-hand (normal-fluid) side are ordered 
from top to bottom with respect to decreasing heat cur- 
rents Q = 160, 130, 100, 70, 40, 20, 0, -20, -40, -60, 
— 70nW/cm^. On the right-hand (superfluid) side they are 
ordered from bottom to top. The horizontal black solid line 
represents the temperature in thermal equilibrium for Q = 0. 
The straight black dashed line represents the superfluid tran- 
sition temperature T\{z). 



and has been applied in Fig. El 

On the left hand side for low altitudes z < — 100/xm 
the system is normal fluid. Here the heat transport equa- 
tion Q = —XtOzT implies that the temperature gradient 
dzT is negative for positive heat currents Q and positive 
for negative heat currents. The values of the gradients 
are considerably large. On the right hand side for high al- 
titudes z > -j-100 fiui the system is superfluid. Here the 
heat is transported convectively following the two-fluid 
model so that the temperatures T(z) are nearly constant 
and the gradients are nearly zero. The intermediate re- 
gion — 100 fiui ^ z -1-100 /im is the superfiuid/normal- 
fiuid interface. Here the temperature profiles interpolate 
the two outer regions. 

For positive heat currents Q > (heat flow upward) 
the slope of the temperature curve T(z) increases without 
a limit on the normal-fluid side for z — >■ —00. However, 
for negative heat currents Q < (heat flow downward) 
the slope increases up to a limiting value which is the 
slope of T\{z) so that in the limit z —go the temper- 
ature profile T{z) approaches a straight line parallel to 
the straight dashed line T\{z). This latter fact is clearly 
observed in the lower left part of Fig. El It represents the 
self-organized critical state predicted by Onuki Q and 
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discovered in the experiment by Moeur et al. 

While Figs. [Hand [21 are calculated for the thermal equi- 
librium at zero heat current Q = 0, we have calculated 
the related curves also for the nonequilibrium state at 
the nonzero heat currents of Fig. [3] Most curves do not 
change very much, the characteristic forms remain qual- 
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itativcly. An exception is the parameter <t defined in Eq. 
p.85p and shown as orange dash-dotted line in Fig. [31 
This parameter is negative in the whole system only for 
small heat currents in the interval — lOnW/cm^ ^ Q ^ 
+20nW/cm^. For larger heat currents outside this inter- 
val the parameter a will change the sign from negative 
to positive at specific altitudes z. For even larger neg- 
ative heat currents Q < — 20nW/cm^ and even larger 
positive heat currents Q > +40nW/cm^ the parameter 
a is positive in the whole system. 



C. Order parameter 

The order parameter in physical units (ip) is calcu- 
lated from the dimensionless renormalized order param- 
eter Y via the renormalization formulas p.l[) and (|3.3ip . 
Putting these equations together and replacing /i — >■ 
we obtain 




(d-2)/2 



(4.4) 



FIG. 4: (Color online) The modulus of the order parameter 
M — as a function of the temperature difference T — 

T\ for the superfluid/normal- fluid interface in gravity. The 
colored solid lines from left to right represent the heat currents 
Q = 160, 130, 100, 70, 40, 20, OnW/cm^ 



Integrating the defining equation p.56p for the zeta func- 
tion C^, we obtain an integral representation for the 
renormalization factor 



(4.5) 



The second equality sign is implied by the fiow-parameter 
transformation (j3.60p together with the running expo- 
nents 1/ and 7], defined in (|3.70p and p.7ip . The upper 
infinite integration boundaries guarantee ~ 1 in the 
limits /i — >■ oo and r — oo which represent the mean- 
field or Gaussian fix point of the RG fiow. If we use the 
correlation length £, — t~'^ we obtain the asymptotic 
formula for the order parameter {tjj) 



with the correct critical exponent /3 defined in p.72p . 

Eqs. (j4.4p and (|4.5p are suited for a numerical calcula- 
tion once the dimensionless renormalized order parame- 
ter Y, the RG fiow parameter r, and the running expo- 
nents (|3.70p - (|3.72p are known. We have calculated the 
order parameter (ip) = M e"'' in the stationary state for 
all those heat currents Q for which we have calculated 
the temperature profiles in the previous subsection. We 
obtain the modulus M and the phase (p of the order pa- 
rameter. Our results for the modulus M are shown in Fig. 
m for positive heat currents Q > (heat fiow upward) . 
The colors of the solid lines correspond to those in Fig. 
[3l Here and in the following figures we omit the lines 
for negative heat currents Q < (heat fiow downward) 
because they make the figures complicated and involved 
but do now show new physics. 

Close to criticality T — T\ the curves are smooth. 
This is an effect of gravity and related to the 
superfluid/normal-fluid interface. The width of the 
smooth region is AT/.g — 25 nK which corresponds to 
the thickness of the interface Azj^g = 200 fim. The 
ratio is approximately the gradient of the superfluid 



transition temperature, i.e. ATj^g/ Azi^g ~ dT\/dz = 
1.273 /iK/cm. For increasing heat currents Q the smooth 
curves are shifted to the left to lower temperatures. This 
fact is related to the depression of the superfiuid tran- 
sition to lower temperatures by a heat current which 
has been observed and investigated in the experiment 
by Duncan, Ahlers, and Steinberg [24 1. 



Away from criticality for lower temperatures T — T\< 
—30 nK the curves approach asymptotically a single line 
which corresponds to the singular order parameter M = 
I ^ {Tx — T)^ for T < T\ in thermal equilibrium 
and zero gravity. In Fig. [3| the asymptotic curves do 
not fall perfectly on a single line. This observation is a 
numerical error in our calculation. In order to stabilize 
the numerical iterations we must add an imaginary part 
to the parameter a defined in p.85p . This imaginary part 
increases with increasing heat current Q and infiuences 
slightly the curves on the superfluid side. 

The physical units cm^^^^ of the order parameter aris- 
ing from the formula (|4.4p for d = 3 dimensions appear 
to be artificial and unphysical. However, since the order 
parameter can not be observed in physical experiments, 
this artifact is not important and no matter of concern. 

The phase of the order parameter ip is dimension- 
less. Its gradient is related to the superfiuid velocity 
Vg — {h/m4)Vip. For nonzero heat currents Q we find 
nontrivial results for the superfluid velocity Vs ■ If we ap- 
proach the interface from the superfluid side, Vg increases 
monotonically. However, on the normal-fluid side, the 
phase If and the superfluid velocity Vg are irrelevant be- 
cause the modulus M approaches zero. 
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FIG. 5: (Color online) The correlation length ^ as a 
function of the temperature difference T — T\ for the 
superfluid/normal- fluid interface in gravity. The colored solid 
lines from left to right represent the heat currents Q — 160, 
130, 100, 70, 40, 20, OnW/cm^. As a reference the black 
dashed line shows the singular correlation length in thermal 
equilibrium and zero gravity. 



D. Correlation length 

The correlation length ^ has been calculated by 
Schloms and Dohm 21[ in thermal equilibrium and zero 
gravity. In the renormalizcd perturbation theory up to 
two-loop order they obtain = M^^J with an am- 
plitude function = 1 + ©(u^). However, since our 
Hartree approximation is first order in u we may ap- 
proximate j4j « 1, so that the correlation length is just 



^ = fi~ ■ This quantity is provided by our numerical cal- 
culation. Our results are shown in Fig. [S] by the colored 
solid lines for the same positive heat currents as in the 
previous figures. In the interface region close to criticality 
T = T\ the colored solid curves are smooth. The corre- 
lation length has a maximum value S,g ~ 50 //m which 
is implied by the gravity acceleration g = 9.81 m/s^ on 
earth. The effect of a small nonzero heat current Q is 
weak. For increasing heat currents Q the position of the 
maximum of the correlation length is shifted slightly to 
lower temperatures. We note that our maximum cor- 
relation length is of the same order of magnitude as 
the characteristic length Jo = 67 yum which was used by 
Ginzburg and Sobyanin [3| within their -i/' theory. 

From Fig. [T] we have inferred the interface thickness 
Azj^g = 200 /im. Thus we calculate the ratio Azj^g/£,g sa 
4 which means that the interface thickness is four times 
the maximum of the correlation length. While in a 
nonequilibrium and/or gravity environment the correla- 
tion length ^ is finite and a smooth function, in equilib- 
rium and zero gravity it shows the well known singular 
behavior ^ ~ jT- TaI"" near criticality for T ^ Tx with 
an exponent v = 0.671. This latter singular correlation 
length is shown by the black dashed line which diverges 



at T = Tx. Far away from criticality which means far 
away from the interface all solid lines converge to a single 
line which is identical with the black dashed line. Thus, 
far away from the interface the gravity g and the heat 
current Q do not have an influence on the correlation 
length ^. Finally, here we do not see an influence of the 
imaginary part of the parameter a we introduce in our 
calculation in order to stabilize the numerical iterations. 



E. Specific heat 

There are two possibilities to calculate the specific 
heat. First, we may calculate the entropy S within 
our renormalization-group theory and then calculate the 
derivative Cx = T{dS/dT)x numerically where any 
quantity X may be ke pt c onstant. This has been done 
in our previous paper [14| where X = Q ot X = VT. 
The entropy S is given by Eqs. (8.10) or (8.12) of Ref. 
[T3- Secondly, we calculate the specific heat directly by 
C' = fcBXo-^mCron wlicrc the renormalizcd specific heat 
Cron is defined in (j3.22p and the renormalization factor 
is defined implicitly in ((3? 



XoZn 



Thus, we obtain 



C 
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( 9p )x. 



} , (4.6) 



a formula which should be compared with the entropy 
(8.10) in Ref. [3 The formula can be simplified if we use 
the asymptotic formulas for the correlation length ^ ~ 
^0 T~'^ and for the coupling parameter 7"^ = (4i'/a)(l — 
6t") where v — 0.671 and a = 2 — dv ^ —0.013 are 
critical exponents and 6 is a known constant. As a result 
we obtain the specific heat 



Cx^B + A{{Aiy/a) + Fx[u]}t- 



together with the amplitude 

Fx[u] 



1 r 

2u 



1 



( 9p ) 



(4.7) 



(4.8) 



This formula should be compared with the entropy (8.12) 
in our previous paper [l3| together with Eqs. (8.13)- 
(8.16). Here A and B are nonuniversal constants which 
can be expressed in the forms A ~ A:B^d/4'fo ^^'^ ^ = 
A^—Aiy/ajb. Alternatively, these constants can be ob- 
tained by fitting the formula to the experimental data 
for liquid *He in a micro gravity environment in space 
[2^ [23} . In this way we obtain A = 2.22 J/molK and 
B = 456 J/molK where the constants are multiplied ad- 
ditionally by the molar volume of liquid ^He at saturated 
vapor pressure [l^] Vx = 27.38 cm'^/mol. 

The amplitude can be compared directly with 

the amplitudes i^±[u] of Dohm Q, if we consider the 
asymptotic limits far away from the interface. The tem- 
perature parameters p and pi are related to each other 
by (|3.36p . In the normal- fluid region far away from the 
interface the renormalizcd order parameter is y = and 
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the amplitudes A and Ai are given by (I3.40p and (|3.4ip . 
The partial derivative can be performed easily so that 
we obtain dp /dpi ~ 1 + 8u Ai ^ 1 — Au which does not 
depend on the variable X that is kept constant. Thus we 
obtain the amplitude 

F+[u] = (2m)-1[1 - 1/(1 - 4u)] = -2 + 0{u) . (4.9) 

In the supcrfluid region pi approaches more rapidly 
than p approaches —1/2. Consequently, in the superfluid 
region far away from the interface the partial derivative is 
dpi /dp = which again does not depend on the variable 
X that is kept constant. Thus we obtain the amplitude 
{2u)~^ . If we compare our results for i^±[M] 
with those of Dohm Q we find agreement for the leading 
terms in powers of u in both cases + and — , respectively. 

We have calculated the specific heat numerically with 
both methods described above using the entropy formula 
(8.12) of our previous paper [ij] and the specific-heat for- 
mula (|4.7p of the present paper. The results agree with 
each other within the accuracy of our Hartree approxi- 
mation which is a self-consistent one-loop approximation 
combined with the renormalization-group theory. This 
agreement is a test for the validity and the accuracy of 
our method presented in this paper. While in the previ- 
ous paper we have used the first method, in this paper 
we prefer the second method, i.e. formula (j4.7|) together 
with (|4.8|) . The reason is that in the present calculation 
the second method provides curves looking smoother and 
more nice. 

Our results are shown in Fig. [6] by the colored solid 
lines for the same heat currents as in the previous figures. 
We have calculated the specific heat Cq where the heat 
current Q and the gravity acceleration g = 9.81 cm/s^ 
are kept constant. Clearly, in the interface region near 
criticality we find smooth curves. The specific heat has 
a maximum slightly below the critical temperature. For 
increasing heat currents Q this maximum is shifted to 
lower temperatures which is related to the depression of 
the superfluid transition temperature observed in the ex- 
periment by Duncan, Ahlers, and Steinberg [IJ]. Fur- 
thermore, the maximum of the specific heat is strongly 
enhanced for increasing heat currents Q. This enhance- 
ment is an effect of the constant heat current Q when 
calculating the specific heat Cn- It has been observed al- 
ready in our previous paper |14l | , where Cq has be calcu- 
lated for the much higher heat current Q = 42.9 ^W/cm^ 
where gravity effects are negligible. The strong enhance- 
ment of the maximum is also compatible with experimen- 
tal measurements of Cq by Harter et al. [2^ . 

Far away from criticality and the interface on both 
sides the colored solid curves converge to a single line, 
respectively. These single lines represent the asymptotic 
specific heat C = B+{A±/a)\t\-"' where t = {T-Tx)/Tx 
is the reduced temperature and a is the critical exponent. 
On the normal-fluid side the single line is perfect. How- 
ever, on the superfluid side it is slightly influenced by the 
imaginary part of the parameter a which we must add 
in our numerical calculation in order to have stable iter- 




""I I I I I I I I — I 

-100 -50 50 100 

T-T^ [nK] 

FIG. 6: (Color online) The specific heat Cq as a function of 
the temperature diflerence T — Tx for the superfluid/normal- 
fluid interface in gravity. The colored solid lines from left 
to right represent the heat currents Q = 160, 130, 100, 70, 
40, 20, OnW/cm^. The black dashed line represents the sin- 
gular specific heat in mi crog ravity fitted to the data of the 
experiment by Lipa et al. [2^ [23| . The black crosses show the 
experimental data for zero heat current in gravity on earth by 
Lipa [2^ where the black dotted line is the related theoretical 
curve for the average specific heat Cq. 



ations. This fact is related to the similar observation in 
our results for the order parameter shown in Fig. SI 

The smooth colored solid lines in Fig. [6] show that the 
critical singularity is rounded by the gravity g and the 
heat current Q. The temperature scale for this round- 
ing is i^Tgj = 25 nK if gravity is the dominating effect. 
We have obtained this value from the thickness of the 
interface Az/^g = 200 ^m. Consequently, the asymptotic 
critical behavior of the speciflc heat and all other sin- 
gular quantities can be observed only for temperatures 
\T ~ T\\ > ATgj = 25 nK away from criticality. Hence, 
the gravity implies that on earth the critical point can 
never be reached. For this reason, experiments to mea- 
sure the asymptotic behavior closer to the critical point 
must be performed in a micro-gravity environment in 
space. 

Lipa et al. [2^ . [23l | have performed a space experiment 
which was called Lambda Point Experiment (LPE) and 
which flew aboard the space shuttle Columbia (STS-52) 
in 1992. They obtained data for the speciflc heat up to 
|r — Ta| = 1 nK. They fitted an asymptotic formula to 
the data and determined the exponent a = —0.013, the 
amplitudes A± and B and some further parameters. The 
resulting fit curve is shown in Fig. [5] by the black dashed 
line. This curve shows the typical lambda of the specific 
heat with a singularity aX T — T\. Away from criticality 
on both sides for \T — T\\ > 50 nK the solid lines and 
the dashed line come close to each other which demon- 
strates the agreement between theory and experiment. 
However, the agreement is not perfect. There remains 
a small discrepancy which is due to the amplitude ra- 
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tio A+ 1 A- because our theory provides an approximate 
value for this amphtude ratio which can never be identi- 
cal to the experimental value. 

We note that in our calculation the specific heat is 
defined locally. It depends on the altitude z so that 
Cq = Cq{z). The rounding of the critical singularity 
is caused by the gradient of T\{z) and by the heat cur- 
rent only. However, in experiments the ^He is in a cell of 
a finite extension. The cell is usually confined by two hor- 
izontal plates at altitudes zi and Z2, where the vertical 
height L — Z2 ~ zi is small and the horizontal exten- 
sions are large. For this reason the experiment measures 
an average specific heat which is defined by the integral 
Cq = {z2 — zi)^^ J^'^ dzCglz). This average process 
smooths the curve additionally. The maximum of the 
average specific heat Cq will be broader than the maxi- 
mum of the related local specific heat Cq . 

In order to minimize the averaging effects the cell 
height L should be chosen as small as possible. How- 
ever, it must be considerably larger than the maximum 
correlation length = 50 because for small L finite 
size effects occur which again smooth and round the crit- 
ical singularity. Consequently, for the cell height L there 
will be an optimum range to obtain best measurements 
on earth. 

In Fig. ini the black crosses represent the experimen- 
tal data of a measurement on earth at zero heat cur- 
rent performed by Lipa [l^. In this case the cell height 
is L = 0.38 mm which causes considerable average ef- 
fects. The maximum of the experimental data is much 
broader than the maximum of the black solid line which 
represents the local specific heat Cq for Q = 0. We 
have calculated the related average specific heat Cq for 
(5 = which is shown by the black dotted line. This 
latter curve shows a much broader maximum at critical- 
ity which agrees with the experimental data. Since the 
ratio L/£^g = 7.6 is large, finite size effects are small. 
However, the Dirichlet boundary conditions of the order 
parameter at the cell walls imply that nevertheless the 
finite size effects cause a depression of the data which is 
clearly observed in Fig. |5] because the experimental data 
(black crosses) are below the theoretical curve (black dot- 
ted line). Thus we conclude that the experimental data 
obtained in a measurement in gravity on earth agree with 
our theory. 



F. Thermal conductivity and resistivity 

The thermal conductivity At is defined locally by the 
heat transport equation Q = — AtVT. Resolving this 
equation we obtain the thermal conductivity explicitly 
as At = |Q|/|Vr|. Next we insert the renormalization 
equations for the heat current (|4.3p and for the temper- 
ature gradient p.65p . Thus, we obtain 



The dimensionless heat current q and the dimensionless 
renormalizcd gradient ^DAp are variables in our numer- 
ical calculation. Hence Eq. (j4.10p is well suited for an 
explicit calculation of the thermal conductivity. 

The dimensionless renormalized heat current is defined 
in p.82[) . Far away from the interface in the normal- 
fluid region the order parameter Y and hence the last 
term of (|3.82p is zero. On the other hand, the amplitude 
Ai in the first term of (|3.82l) reduces to Ai = —1/2 
following (|3.41D . Thus Eq. (I4.10p provides the thermal 
conductivity 
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in the normal-fluid region. This result is well known and 
agrees with the linear-response calculation of Dohm 0] . 
Since our result is derived in Hartree approximation, the 
agreement is up to one-loop order. 

Far away from the interface in the supcrfluid region the 
heat current Q is nonzero where the temperature gradi- 
ent VT is zero. Hence the thermal conductivity At is 
infinite. This fact is also seen in the formulas (j3.82|) and 
(|4.10p . A zero dimensionless renormalized temperature 
gradient ^DAp = implies a zero first term in p.82p and 
a zero denominator in (j4.10[) . On the other hand, the sec- 
ond term in (j3.82[) is minus the dimensionless renormal- 
izcd supcrfluid current which is nonzero. Consequently, 
Eq. (j4.10|) provides an inflnite thermal conductivity once 
again in the supcrfluid region. 

A related quantity is the thermal resistivity pT = 1/At 
which is the inverse of the thermal conductivity. In the 
normal-fluid region the thermal resistivity pT is finite. 
On the other hand, in the superfiuid region it is zero. 
For this reason the thermal resistivity pT is well suited 
for a graphical representation. In Fig. [7]the results of our 
numerical calculation are shown as colored solid lines for 
several positive values of the heat current Q. These heat 
currents are the same as those in the previous figures. 
Clearly, in the interface region near criticality the colored 
solid curves are smooth lines. The critical singularity at 
T = Ta is smoothed by the nonzero values of gravity g 
and the heat current Q. For increasing heat currents Q 
the colored solid lines are shifted to lower temperatures. 
This fact is again related to the depression of the su- 
pcrfluid transition temperature by nonzero heat currents 
following Duncan, Ahlers, and Steinberg [2^ . 

On both sides far away from criticality T = T\ and 
from the interface the colored curves asymptotically ap- 
proach single lines, respectively. In the normal-fluid re- 
gion the asymptotic thermal resistivity is given by (|4.1ip , 
where in the supcrfluid region it is just zero. 

The local thermal resistivity pT has been measured 
for the supcrfluid/normal-fluid interface in the exper- 
iment by Day et al. [23|. While gravity on earth is 
g = 9.81 m/s^, the heat current is flowing upward from 
bottom to top so that Q is positive. The values of heat 
current Q are the same as in our calculation. For this 
reason, our colored solid lines can be compared directly 
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FIG. 7: (Color online) The thermal resistivity as a function of 
the temperature difference T — T\ for the superfluid/normal- 
fiuid interface in gravity. The colored solid lines from left to 
right represent the heat currents Q = 160, 130, 100, 70, 40, 
20, OnW/cm^. The data points represent the experimental 
data by Day et al. [23] taken for the same heat currents except 
Q = 0. Related to the theoretical curves the data are ordered 
from left to right. 

with the experimental data. In Fig. [7] the experimental 
data for several heat currents Q are represented by points 
of several symbol types. The data points can be related 
to the solid curves by their color or alternatively by there 
order from left to right. An exception is the black solid 
line for zero heat current Q = 0. In this latter case no 
experimental data are available. However, the black solid 
line is very close to the red solid line for Q = 20 nW/cm-^. 
Consequently for Q < 20nW/cm^ the effect of the heat 
current is very small. 

Qualitatively, the experimental data agree with the 
theoretical solid lines. However, there are quantita- 
tive discrepancies. First of all, in the normal-fluid re- 
gion well above criticality for temperature differences 
T — Tx > 50 nK the thermal resistivity converges to a 
single line represented by the formula (|4.1ip . We have 
multiplied the theoretical results with a correction factor 
which slightly differs from unity. In this way we achieve 
that far away from criticality the experimental data are 
lying on top of the theoretical curves. This correction 
factor is justified, because in Ref. the model- param- 
eters were adjusted for the theoretical specific heat and 
thermal resistivity in two-loop order where in the present 
paper and also in our previous papers [isl [l^ the quanti- 
ties were calculated in the Hartree approximation which 
is a self-consistent one-loop approximation. 

For larger heat currents Q and temperatures T slightly 
below Tx our solid lines show some bumps which are 
probably artifacts of our approximation. The magnitude 
of the artifacts is within the accuracy of our approach. 
In the superfluid region well below criticality for temper- 
ature differences T — Ta < — 50nK the thermal resistiv- 
ity pt is very close to zero both in theory and experi- 



ment. This fact represents the frictionlcss heat transport 
by the superfluid/normal- fluid counterflow in superfluid 
^He. Vortices are not present in our calculation presented 
here. 



G. Second sound and time dependent phenomena 

Until now we have considered only stationary nonequi- 
librium states where a constant heat current Q = Q 
is flowing vertically in the ''He and where all time- 
dependent phenomena are relaxed. However, our numer- 
ical calculation solves the time-dependent model- equa- 
tions so that time-dependent phenomena can be treated 
explicitly. In this case the order parameter {Tp{z,t)), 
the temperature T(z,t), and the heat current density 
Q(z,t) = Q{z,t)ez are functions of altitude z and time 
t. In superfluid ''He the most important dynamic phe- 
nomenon is second sound. We can generate a second 
sound pulse on the upper cell boundary Z2 if in the ex- 
ternal heat source (j4.ip we choose a time-dependent up- 
per source function Q2{t). Then the phase of the order 
parameter ip(z,t), the temperature T{z,t), and the local 
heat current Q{z,t) show a pulse which is traveling down- 
ward toward the superfluid/normal- fluid interface. Ap- 
proaching gradually the interface the width of the pulse 
increases. Once the interface is reached the pulse disap- 
pears by broadening where nearly nothing is reflected. In 
the end the second-sound pulse is absorbed nearly com- 
pletely by the interface. 

A heat pulse can be generated also on the lower cell 
boundary zi if the lower source function Qi(t) is chosen 
time dependent. However, in this case the response of 
the system is less spectacular because on the normal-fluid 
side the heat is transported diffusively. Nevertheless, our 
time-dependent investigations yield an important result. 
All time-dependent perturbations of the stationary states 
with a constant heat flow relax and disappear. Thus 
we conclude that the superfluid/normal-fluid interface in 
gravity and in the presence of a vertical heat flow is a 
stable physical configuration. 

V. COMPARISON WITH OUR PREVIOUS 
APPROACH FOR LARGE HEAT CURRENTS 

A. Solutions 

We have found two different solutions of the model-F 
equations for superfluid ''He in the nonequilibrium state 
where a heat current is flowing. In the first case the order 
parameter (ip) = M e*"^ is nonzero and no vortices are 
present. The heat is transported convectivcly without 
any friction by the superfluid/normal-fluid counterflow 
so that the temperature gradient is zero. This solution is 
investigated in the present paper. In the second case the 
order parameter is zero due to fluctuations of the phase 
(p by moving vortices and quantum turbulence. Here the 
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moving vortices imply a small thermal resistivity which 
causes a small temperature gradient. This latter solution 
was investigated in our previous papers [H, [11] . The two 
different solutions exist only in the superfluid and inter- 
face region where a nonzero order parameter is possible. 
In the normal-fluid region the solution is unique because 
here the order parameter is always zero. 

The two solutions are controlled by the parameter a 
defined in Eq. p.85|) . The first solution exists for nega- 
tive and positive values of cr where the second solution 
exists only for positive a. For the superfluid/normal-fluid 
interface in gravity on earth this means that the first so- 
lution exist for small and large heat currents. On the 
other hand the second solution exists only for large heat 
currents |Q| > 70nW/cm^ where the heat current is the 
major and gravity is the minor influence. However, the 
numerical calculations of our present paper show that in 
practice the first solution is stable only for small heat cur- 
rents in the interval — 70nW/cm^ ^ -1-160 nW/cm^. 
For larger heat currents outside the interval we do not 
obtain a stable solution. We do not know if this instabil- 
ity is a property of our iteration procedure only which we 
describe in Subsec. IIII Fl However we guess that there is 
a physical instability beyond a certain critical heat cur- 
rent. This means we expect a discontinuous first-order 
like transition in the noncquilibrium state near the heat 
current Qq w 70nW/cm^ which separates the gravity 
dominated regime (first solution) from the heat current 
dominated regime (second solution). Nevertheless, there 
will be an overlap region where both solutions exist. 



B. Stability 

The time-dependent nature of our numerical calcula- 
tion provides a test for the physical stability of our sta- 
tionary solutions. We may add a small perturbation to 
the solution and then start the calculation. After a time 
difference St 2 s the perturbations relax and disappear 
so that the system returns to the stationary state. We 
find this behavior for all stationary solutions with a con- 
stant heat flow which we have presented in this paper. 
Thus, we conclude that the first-type solutions are sta- 
ble which describe physical states with a nonzero order 
parameter and no vortices. This result is expected be- 
cause the order-parameter variation and second sound 
are damped. We note that we have a numerical instabil- 
ity for larger heat currents. However, this latter insta- 
bility is unphysical and has a completely different nature 
because it appears on the short length scales Az of the 
discretization of the altitude coordinate. 

In our previous papers p^ . we did not prove the 
stability of the second-type solutions which describe su- 
perfluid states with vortices and a zero average order pa- 
rameter. The reason is that in our previous papers we 
did not solve the model-F equations as partial differen- 
tial equations. Consequently, we do this now and pro- 
vide the proof in the following. However, we note that in 



this case the stability is nontrivial because the system is 
in a superfluid state where the average order parameter 
is (ip) ~ and the temperature is T < T\. We solve 
the time-dependent renormalized model- equations for 
a second-type solution with a small perturbation. For 
simplicity we consider a self-organized critical state with 
a constant heat Q current and constant temperature gra- 
dients VT = VTa because this state is spatially homoge- 
neous before the perturbation is applied. Unfortunately, 
the state is not periodic because we cannot require peri- 
odic boundary conditions for all quantities. Exceptions 
are the temperatures T[z.t), T\{z) and the phase of the 
order parameter (p{z,t). Nevertheless, we can general- 
ize the boundary conditions. For the latter three quan- 
tities we require pseudoperiodic boundary conditions in 
the sense of the impossible objects of the famous Dutch 
graphic artist M.C. Escher [g^, i.e. 

T{z + L,t) = Tiz,t) + ATL , (5.1) 
Tx{z + L) = Tx{z) + ATl , (5.2) 
ip{z + L,t) = ip{z,t) + AipL{t) . (5.3) 

Since the critical temperature T\(z) is linear in z, from 
Eq. we obtain ATl = {d^Tx)L. On the other 

hand from Eq. (jS.ip we obtain ATl = ^{Q/^t)L where 
At is the thermal conductivity. We note that for the 
self-organized critical state both results for the period 
constant ATl must be equal. This fact implies Q = 
—Xt^OzTx). The period constant of the order-parameter 
phase AifLit) can be determined by investigating the 
renormalized model-F equation (|3.79|1 . On the right- 
hand side we replace Ap{z, t) in favor of the temperature 
T{z,t) by inserting (|3.27|1 . We insert the complex order 
parameter Y = r] e"^ and derive an equation for the or- 
der parameter phase ip{z,t). We consider this equation 
for the altitudes z = zq and 2 = 2:0 + ^ and then sub- 
tract the resulting equations. Thus, as a result we obtain 
dt AipL{t) = go{ATL/T\). We integrate this equation and 
obtain the period constant 

A^L{t) - A^L{to) + gQ{ATL/Tx){t - to) . (5.4) 

Clearly, the period of the order-parameter phase depends 
linearly on the time t. We conclude and find that for the 
self-organized critical state with small perturbations the 
renormalized model- _F equations can be solved numeri- 
cally using Escher pseudoperiodic boundary conditions 
which are defined by Eqs. (|5.ip - (|5.3p . Our numerical test 
provides the following result. Any small perturbation of 
the self-organized critical state relaxes and disappears af- 
ter a time difference of (5i w 2 s. Thus we conclude that 
the second-type solutions are stable. This means that the 
noncquilibrium states considered in the previous paper 
arc stable. 

A matter of special interest is the relaxation of the 
order parameter. We start at time to with a small con- 
stant perturbation Y{z,to) = 7706*'^°. This means 770 is 
small but nonzero, and (po is constant. We consider the 
renormalized model--F' equation (j3.79p which describes 
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the time evolution of the order parameter Y{z,t) = 
77(2, e*'^*^^'*^. Decomposing the equation with respect 
to the modulus ri{z,t) and the phase ip{z,t) from the 
first term on the right hand side we infer the damping 
for the modulus 

D = g„{2-fT){w'/F)[p, + i^V^f] . (5.5) 

This damping is an inverse relaxation time. The solution 
of the model- _F equation is stable whenever this damp- 
ing is positive and unstable otherwise. The prefactors arc 
always positive so that the crucial quantity is the expres- 
sion in the square brackets. For times t shortly after the 
beginning of the calculation to the phase is expected to 
be ip{z,t) « (po so that T/ip « 0. Consequently, the main 
contribution is the dimensionless modified temperature 
parameter pi. The stability of the solution depends on 
its sign. In Fig. [T] pi is plotted as the blue dashed line 
for the superfluid/normal- fluid interface. A similar curve 
is obtained for the self-organized critical state if pi is 
plotted as a function of the heat current Q. While in the 
normal-fluid region pi is positive, in the superfluid region 
it is negative. For the heat current Q = 170nW/cm-^ we 
find the minimum value pi^min = —0.39. Thus, we con- 
clude that the solution is stable in the normal-fluid region 
but unstable in the superfluid region. 

However, the instability is true only for short times 
where t — to is small. For longer times we must investigate 
the space and time dependence of the order-parameter 
phase ip{z,t). From Eq. (|5.4p and from the renormalizcd 
model- equation p.79p we infer 

ipiz,t) = (fio + AipLit) {z ~ zo)/L 

= ifo + go {ATL/Tx){t- to) {z ^ zo)/ L . (5.6) 

Here zo is the altitude where the temperature equals the 
reference temperature, i.e. T{zo,to) = Tq. For the last 
equality sign we have used (j5.4p realizing that in our 
case at the beginning t ^ to the Escher period constant 
is A(pL{to) = 0. Differentiating with respect to the al- 
titude coordinate z and multiplying by ^ we obtain the 
dimensionless gradient of the order-parameter phase 

^Vcp = go [ATl/Tx) [t - to) {i/L) e, . (5.7) 

Consequently, we find (^V(p)^ ^ [t — to)^ ■ This means 
that in Eq. (|5.5p the second term in the square bracket 
increases with time. Even though at the beginning the 
square bracket may be negative because of the first term, 
after a short time the second term makes the square 
bracket positive. Thus we conclude: Even though there 
may be an instability for short times, the increase of the 
gradient of the phase (j5.7p makes the damping (j5.5p fi- 
nally positive so that the time evolution of the order pa- 
rameter is flnally stable. 

The stability of the second model- equation p.80p is 
easily proven. We insert the dimensionless renormalizcd 
heat current (j3.82|) and neglect its last term because it 
is squared in the small nonzero order parameter. As a 



result we obtain a diffusion equation for the dimension- 
less renormalizcd temperature parameter Ap{z,t). Since 
the related diffusion constant is positive, this equation is 
always stable. The covariant derivatives in this equation 
do not affect the stability. 

We summarize that we have presented an explicit proof 
for the stability of the self-organized critical state which 
is a spatially homogeneous second-type solution of the 
model-i^ equations. We expect that also the more gen- 
eral second-type solutions for spatially inhomogeneous 
systems arc stable, which describe the superfluid / normal- 
fluid interface in our previous papers Our numer- 

ical calculations of the present paper support this expec- 
tation. However, we note that the stability is nontrivial 
because the system is in a superfluid state where the av- 
erage order parameter is = and the temperature is 
T <Tx. 



VI. DISCUSSION AND CONCLUSION 

Onuki 0, Q and and later Weichman and Miller ^ 
have also investigated the superfluid/normal-fluid inter- 
face within model F . They obtain temperature profiles 
which agree qualitatively with our results shown in Fig. 
O However, they did not use the rcnormalization-group 
theory and the related coupling parameters which have 
been determined by Dohm [8|. For this reason, it is not 
possible to compare the results quantitatively. Weichman 
and Miller furthermore considered the self-organized 
critical state for a heat current flowing downward. While 
in the superfluid region usually the temperature profile 
is flat, they obtain phase slips in the order parameter 
which produce a stair-case like temperature profllc. In 
this way they obtain a temperature gradient VT which 
on average equals the gradient VT^ as required for the 
self-or gan ized critical state. More recently Yabunaka and 
Onuki [30| performed a three-dimensional numerical sim- 
ulation based on model F in order to investigate the self- 
organized critical state and the superfluid/normal-fluid 
interface. They observed the formation and motion of 
vortices and phase slips which produce a nonzero tem- 
perature gradient VT on average in the superfluid region 
which compensates VTx. 

A sophisticated theory for mutual friction, quantum 
turbulence, and the dynamics of vortices in superfluid 
■*He was developed long time ago by Vinen jsj. A mea- 
sure for the quantum turbulence is the density of the 
vortices which is deflned as the total length of the vortex 
lines per volume. For this vortex density a rate equation 
is derived. On the right hand side of this equation there 
is a term for the generation and a term for the decay of 
vortices and quantum turbulence. Vortices are usually 
generated by a nucleation process. This means that an 
energy barrier must be overcome which strongly reduces 
the generation rate. 

We believe that model-T includes the effects of vor- 
tices and quantum turbulence correctly so that the Vi- 
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nen theory can be derived if the model-F equations can 
be solved exactly without any approximation. However, 
our two solutions which are derived within the Hartree 
approximation are idealized solutions of model F. They 
describe the two phases of a first-order nonequilibrium 
transition but do not include metastability and the nu- 
cleation process. 

For the check between theory and experiment an im- 
portant quantity is the thermal resistivity px in super- 
fluid ^He for T < T\ and larger heat currents Q > 
lOOnW/cm^ induced b y th e effect of vortices. In Fig. 
3 of our previous paper |14l | we have compared the result 
of our second-type solution with experimental data by 
Baddar et al. [321 • The experimental thermal resistivity 
is lower by a factor of 20 than our theoretical result. A 
plausible explanation for this discrepancy is the follow- 
ing. In the experiment the heat current is flowing upward 
from bottom to top. The related superfluid current is 
flowing in the opposite direction, i.e. downward from top 
to bottom. Consequently, the vortices are transported 
together with the superfluid current downward. Since 
the order parameter (ip) increases and the vortex density 
decreases with the altitude z, the downward transport of 
the vortices together with the metastability of the nu- 
cleation process can reduce the vortex density consider- 
ably. Since the hopping over energy barriers in the nu- 
cleation process causes exponential factors, a reduction 
of the thermal resistivity by a factor 20 is plausible. 

The situation is different for the self-organized criti- 
cal state. Here the system is spatially homogeneous be- 
cause the heat current Q and the temperature difference 
T(z) — T\{z) = AT are constant and do not depend the 
altitude z. The nucleation process reaches an equilibrium 
state, so that our theory will predict the density of the 
vortices and the related thermal resistivity correctly. A 
flrst agreement with the theory was found in the experi- 
ment by Moeur et al. [l^. Here the temperature differ- 
ence AT was measured as a function of the heat current 
Q. The comparison with our theory is shown in Fig. 3 of 



our previous paper [T^. The agreement is encouraging. 
However, for very large heat currents Q > 2/iW/cm^ 
deep in the superfluid region a deviation was found. This 
deviation may be a problem of the temperature measure- 
ment because the temperature T is never measured in the 
bulk of the system but always on the surface. In order 
to avoid the Kapitza resistance which implies a temper- 
ature jump on the surface, the temperature is usually 
measured by thermometers on the side walls. However, 
in cases where vortices are present, a superfluid flow par- 
allel along a side wall may cause a transverse Kapitza 
resistance so that there is a temperature jump also on a 
side wall. 

In a recent experiment Chatto et al. per- 
formed an experiment to measure the thermal conductiv- 
ity/resistivity indirectly where the explicit measurement 
of the temperature is avoided. Instead they measured 
the velocity v oi a propagating thermal mode as a func- 
tion of the heat current Q in the interval 30nW/cm^ < 
Q ^5 15/iW/cm^. On the other hand, they derived a 
theoretical curve for the velocity v from our theoretical 
results for the thermal conductivity [l^, [ij] . They find 
very good agreement between theory and experiment in 
the whole range of heat currents Q, even for the largest 
values Q ~ 15 /iW/cm^ deep in the superfluid region. 
We conclude that this experiment is an important ver- 
ification of our theory. This means our theory [H Q 
describes the effects of vortices, mutual friction, and the 
thermal conductivity /resistivity correctly on a quantita- 
tive level for the self-organized critical state. 



Acknowledgments 

We would like to thank Prof. Dr. V. Dohm for helpful 
comments on the manuscript and Prof. Dr. J. A. Lipa for 
experimental data of the specific heat. 



[1] G. Ahlers, Phys. Rev. 171, 275 (1968); J. Low Temp. 

Phys. 84, 173 (1991). 
[2] V.L. Ginzburg and A.A. Sobyanin, Sov. Phys. Usp. 19, 

773 (1976). 

[3] A. Onuki, J. Low Temp. Phys. 50, 433 (1983); J. Low 

Temp. Phys. 55, 309 (1984). 
[4] A. Onuki Jpn. J. Appl. Phys. 26, 365 (1987). 
[5] P.O. Hohenberg and B.I. Halperin, Rev. Mod. Phys. 49, 

435 (1977). 

[6] B. L Halperin, P.C. Hohenberg, E.D. Siggia, Phys. Rev. 

Lett. 32, 1289 (1974); Phys. Rev. B 13, 1299 (1976); 

Erratum Phys. Rev. B 21, 2044 (1980). 
[7] V. Dohm, Z. Phys. B 60, 61 (1985); Z. Phys. B 61, 193 

(1985). 

[8] V. Dohm, Phys. Rev. B 44, 2697 (1991); Phys. Rev. B 

73, 092503 (2006); Phys. Rev. B 73, 099901(E) (2006). 
[9] A. Singsaas and G. Ahlers, Phys. Rev. B 29, 4951 (1984); 



Phys. Rev. B 30, 5103 (1984). 
[10] W. Y. Tarn and G. Ahlers, Phys. Rev. B 32, 5932 (1985); 

Pliys. Rev. B 33, 183 (1986). 
[11] R. Haussmann and V. Dohm, Phys. Rev. Lett. 67, 3404 

(1991) ; Z. Phys. B 87, 229 (1992). 

[12] R. Haussmann and V. Dohm, Phys. Rev. B 46, 6361 

(1992) . 

[13] R. Haussmann, J. Low Temp. Phys. 114, 1 (1999). 

[14] R. Haussmann, Phys. Rev. B 60, 12349 (1999). 

[15] W.A. Moeur, P.K. Day, F.-C. Liu, S.T.P. Boyd, M.J. 

Adriaans, and R.V. Duncan, Phys. Rev. Lett. 78, 2421 

(1997). 

[16] R. Haussmann, unpublished notes (2008); A. Nikodem, 
Diploma Thesis, University of Konstanz, Konstanz, Ger- 
many (2008). 

[17] A.R. Chatto, R.A.M. Lee, R.V. Duncan, and D.L. Good- 
stein, J. Low Temp. Phys. 148, 519 (2007). 



24 



[18] L.P. Pitaevski, Sov. Phys. JETP 8, 282 (1959). 

[19] A.L. Fetter and J.D. Walecka, Quantum theory of many- 
particle systems (McGraw-Hill, New York 1971). 

[20] R. Haussmann, Phys. Rev. B 60, 12373 (1999). 

[21] R. Schloms and V. Dohm, Nucl. Phys. B 328, 639 (1989). 

[22] J. A. Lipa, D.R. Swanson, J. A. Nissen, T.C.P. Chui, and 
U.E. Israelsson, Phys. Rev. Lett. 76, 944 (1996). 

[23] J. A. Lipa, J. A. Nissen, D.A. Strieker, D.R. Swanson, and 
T.C.P. Chui, Phys. Rev. B 68, 174518 (2003). 

[24] R.V. Duncan, G. Ahlers, and V. Steinberg, Phys. Rev. 
Lett. 60, 1522 (1988). 

[25] A.W. Harter, R.A.M. Lee, A. Chatto, X. Wu, T.C.P. 
Chui, and D.L. Goodstein, Phys. Rev. Lett. 84, 2195 
(2000). 

[26] J. A. Lipa, private communication (2011). 



[27] P.K. Day, W.A. Moeur, S.S. McCready, D.A. Sergatskov, 
F.-C. Liu, and R.V. Duncan, Phys. Rev. Lett. 81, 2474 
(1998). 

[28] D.R. Hofstadter, Godel, Escher, Bach: An Eternal 
Golden Braid (Basic Books, New York 1979). 

[29] P.B. Weichman and J. Miller, J. Low Temp. Phys. 119, 
155 (2000). 

[30] S. Yabunaka and A. Onuki, Phys. Rev. B 82, 024501 
(2010). 

[31] W.F. Vinen, Proc. R. Soc. A 240, 114 and 128 (1957); 

242, 493 (1957); 243, 400 (1958). 
[32] H. Baddar, G. Ahlers, K. Kuehn, and H. Fu, J. Low 

Temp. Phys. 119, 1 (2000). 



